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Abstract 

Accurate  and  efficient  tire  models  for  soft  soil  operations  are  needed  for  mobility  simulations  of 
wheeled  vehicles  operating  in  off-road  conditions.  The  tire  model  is  essential  in  a  vehicle  model; 
the  accuracy  of  the  predicted  tire  forces  and  moments  has  a  large  impact  on  assessing  the  vehicle 
performance,  reducing  the  number  of  stuck  vehicles,  avoiding  rollovers,  and  developing  advanced 
control  strategies.  While  finite  element  method  leads  to  the  most  detailed  tire-soil  interaction  mod¬ 
els,  their  complexity  and  extensive  computational  effort  make  them  less  than  ideal  for  the  applica¬ 
tions  envisioned,  requiring  hundreds  of  thousands  of  evaluations.  Semi-analytical  models  such  as 
the  F-tire  were  developed  for  on-road  applications.  The  semi-analytical  soft-soil  tire  model  devel¬ 
oped  in  this  study  employs  a  similar  construction  as  F-tire  and  adds  an  advanced  tire-soil  contact 
model.  The  tire  is  discretized  in  a  user-defined  number  of  lumped  masses  structured  in  three  layers 
(two  for  the  sidewall,  one  for  the  tread  and  belt),  connected  by  elastic  and  damping  elements  con¬ 
nected  in  various  configurations.  Simulation  results  for  several  case  studies  performed  with  the  tire 
model  developed  in  this  study  are  presented. 

In  this  paper  we  also  present  the  experimental  work  related  to  soil  testing  and  tire  instrumenta¬ 
tion  and  testing,  as  needed  to  collect  input  data  and  to  validate  the  new  off-road  tire  model.  The  soil 
of  interest  is  silty  sand.  Soil  characterization,  performed  in  collaboration  with  Schnabel  Engineer¬ 
ing,  will  be  presented.  A  physical  tire  has  been  instrumented  with  sensors  able  to  collect  deflection 
information  in  real  time.  Such  information,  in  addition  to  the  information  received  from  a  string 
potentiometer,  is  critical  in  allowing  the  estimation  of  the  true  sinkage  (versus  just  measuring  the 
tire  rut).  Experimental  results  obtained  for  tire  deflection  with  the  instrumented  tire  on  a  quarter-car 
test  rig  will  also  be  presented.  Furthermore,  since  ensuring  controlled,  repeatable  testing  conditions 
is  a  very  important  aspect  for  data  collection,  the  efforts  made  to  guarantee  that  the  normal  load 
applied  to  the  tire  will  stay  constant  during  testing  will  also  be  discussed.  The  influence  of  tire  and 
vehicle  parameters  and  soil  characteristics  on  the  tire  dynamics  is  also  investigated. 
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1  Model  Overview 

The  proposed  soft  soil  tire  model  is  in¬ 
tended  for  mobility,  traction,  and  handling 
applications,  and  it  can  also  be  extended  for 
durability  studies.  For  such  reasons,  a  com¬ 
plex  semi-analytical  tire  modelling  approach 
is  the  best  fit. 


The  schematic  of  the  proposed  model  can 
be  observed  in  Fig.  1.  It  employs  a  discretized 
lumped  mass  approach  that  uses  springs  and 
dampers  in  multiple  configurations  to  repre¬ 
sent  the  different  sections  of  the  tire.  One 
novelty  of  this  model  is  that  it  is  structured  in 
three  parallel  planes,  two  of  them  representing 
the  sidewalls  (pictured  in  blue),  and  one  of 
them  represents  the  belt  and  tread  (pictured  in 
red)  of  the  tire.  By  differentiating  between  the 
sidewall  and  the  belt  a  more  realistic  applica- 
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tion  of  the  local  forces  can  be  implemented,  as 
can  be  seen  in  Fig.  1.  Each  plane  consists  of 
Nm  number  of  masses;  Nm  is  suggested  be¬ 
tween  80  and  100,  but  it  is  user-defined  in  an 
effort  to  increase  the  versatility  of  the  model. 
This  provides  the  user  with  flexibility  in  de¬ 
termining  the  resolution  of  the  model. 

The  rigid  wheel  represents  the  rim  of 
the  tire  and  is  directly  connected  to  the  side- 
wall  layers.  Thus,  the  wheel  has  six  degrees  of 
freedom;  three  translational  and  three  rota¬ 
tional.  On  the  other  hand,  each  lumped  mass 
has  three  degrees  of  motion:  the  translational 
motion  in  all  directions.  Moreover,  there  is 
relative  motion  between  masses  in  the  same 
plane  and  between  masses  in  adjacent  planes, 
as  well.  Relative  motion  is  also  allowed  be¬ 
tween  lumped  masses  and  the  rigid  wheel  in 
the  circumferential  direction,  which  is  another 
novelty  introduced  by  this  model.  The  total 
number  of  degrees  of  freedom  of  the  model 
are, 


Tycp,s  T yCp,bt 


DOF  = 


1.1  Sidewall  Element 

The  sidewall  element  is  more  complex 
than  the  belt  and  tread  element  because  it  is 
directly  connected  to  the  rigid  wheel.  A  dia¬ 
gram  of  the  model  can  be  observed  in  Fig.  2 
and  Fig.  3,  and  a  description  of  each  element 
can  be  found  in  Table  1.  Notice  how  each 
sidewall  mass  has  elastic  and  damping  ele¬ 
ments  in  all  three  directions.  The  mass  of  a 
single  sidewall  element  (mis)  is  calculated 
using  Eq.  (2). 

Q-5msidewall 

m,s  =  — jj -  m 

iV77l 

The  undeformed  radius  of  each  side- 
wall  element  is  the  undeformed  radius  of  the 
tire  minus  the  tread  height  as  shown  below, 

Rs  =  Ru-h  (3) 


Fig.  2:  Sidewall  diagram.  Left-  the  in-plane  connec¬ 
tions;  right  -  the  out-of-plane  connections.  (k8  and  c6  not 
pictured ) 


Fig.  1:  Diagram  of  the  tire  model.  Note  that  all  the  soil 
forces  are  labeled  in  black. 

The  tire  model  requires  several  input 
parameters.  The  soil  model  also  requires  mul¬ 
tiple  input  parameters.  All  the  parameters 
needed  for  the  tire  and  the  soil  model  are 
listed  in  Section  3.1.  The  simulation  inputs 
are:  slip  angle,  slip  ratio,  camber  angle,  tire 
load,  time  for  simulation,  and  finally  either  a 
driving  torque  or  a  longitudinal  velocity.  On 
the  other  hand,  the  outputs  predicted  by  the 
model  are:  tire  deformation,  forces  and  mo¬ 
ments  in  all  directions,  sinkage  for  sidewall 
and  tread  and  belt  layers  and  pressure  distri¬ 
bution  at  the  contact  patch. 


Fig.  3:  Top  view  of  the  elements  in  the  circumferential 
plane. 


1.2  Belt  and  Tread  Element 

The  belt  and  tread  elements  are  different 
from  the  sidewall  element.  This  is  due  to  the 
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fact  that  the  tread  and  belt  layer  is  not  directly 
connected  to  the  rigid  wheel.  The  diagram  for 
the  belt  and  tread  element  can  be  found  in  Fig. 
4.  It  is  important  to  note  that  there  is  an  extra 
set  of  springs  and  dampers  that  connect  the 
sidewall  and  belt  and  tread  layers  in  the  cir¬ 
cumferential  direction  that  is  not  shown  in  this 
figure.  However,  the  extra  set  can  be  seen  in 
Fig.  3.  A  description  of  each  element  can  be 
found  in  Table  1. 


Fig.  4:  Belt  and  tread  element  diagram.  Left  -  in-plane 
view;  right  -  out-of-plane  view.  (k8  and  c8  not  pictured) 

It  is  also  important  to  note  that  the  un¬ 
deformed  radius  of  the  belt  and  tread  plane 
(. Ru )  is  equal  to  what  is  traditionally  under¬ 
stood  by  the  undeformed  radius  of  the  tire 
(. Ru );  this  is  larger  than  the  undeformed  radius 
of  the  sidewall  layer;  thus  the  model  accounts 
for  the  curved  shape  of  the  tire. 


Table  1,  Description  of  tire  model  elements 


Element 

Definition 

ki  &  ci 

Sidewall  radial  spring  stiffness 
and  damping  coefficient  (in¬ 
plane) 

k2  &  c2 

Wheel-sidewall  circumferential 
spring  stiffness  and  damping 
coefficient  (in-plane) 

k3 

Inter-element  radial  spring  stiff¬ 
ness  (in-plane) 

ki  &  C4 

Circumferential  inter-element 

spring  stiffness  and  damping 
coefficient  (in-plane) 

k5  &  c3 

Lateral  inter-element  spring 
stiffness  and  damping  coefficient 
(out-of-plane) 

k6  &  c5 

Radial  inter-element  spring  stiff¬ 
ness  and  damping  coefficient 
(out-of-plane) 

kp  &  C7 

Sidewall-wheel  lateral  spring 
stiffness  and  damping  coefficient 
(out-of-plane) 

k8  &  c6 

Circumferential  inter-element 

spring  stiffness  and  damping 
coefficient  (out-of-plane) 

The  mass  of  a  single  belt  and  tread  element 
(; mi,bt )  is  calculated  using  the  equation, 

^ tread  5"  Albeit  (4) 


1.3  Coordinate  systems 

The  model  uses  the  following  coordi¬ 
nate  systems  to  describe  the  position  of  the 
wheel  and  of  the  lumped  masses: 

The  X-Y-Z  coordinate  frame  is  the  global  ref¬ 
erence  frame.  This  global  coordinate  system 
follows  the  ISO  8855  standard,  presented  in 
Fig.  5.  Therefore,  the  positive  z-axis  is  point¬ 
ing  upwards. 


Fig.  5:  ISO  8855  wheel  coordinate  system.  Adapted 
from  [1]. 

The  xo-yo-Zo  is  the  local  reference  frame  fixed 
to  the  rigid  wheel. 
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Each  lumped  mass  can  be  positioned 
with  respect  to  the  local  wheel  reference 
frame.  Since  the  lumped  masses  have  relative 
displacements  in  all  three  directions  with  re¬ 
spect  to  the  rigid  wheel,  the  springs  and 
dampers  included  in  the  model  along  those 
directions  will  sustain  deformations/ 
displacement.  To  ease  the  calculation  of  de¬ 
formation  forces  in  the  springs  and  damping 
forces  in  the  dampers,  these  forces  are  com¬ 
puted  along  the  directions  of  a  local  reference 
frame,  Xi-yt-zt ,  attached  to  each  of  the  lumped 
mass,  as  seen  in  Fig.  6.  The  lumped  mass 
frames  are  oriented  in  such  a  way  that  the  x, 
reference  frame  points  in  the  circumferential 
direction  and  the  Zi  reference  points  in  the 
radial  direction.  Thus,  it  is  worth  noticing  that 
these  lumped  mass  reference  frames  are  in 
fact  always  rotated  with  a  fixed  angle  with 
respect  to  the  directions  of  the  position  vec¬ 
tors  that  locate  the  mass  in  the  rigid  wheel 
reference  frame. 

The  simulation  is  carried  out  in  two 
steps,  as  explained  below.  To  illustrate  the  use 
of  the  coordinate  systems,  two  random  masses 
where  chosen  in  a  two-dimensional  space,  the 
first  one  is  denoted  with  a  letter  C  and  it  is 
mass  1.  The  second  is  mass  i  and  it  is  denoted 
by  the  letter  A.  Initially,  the  masses  are  locat¬ 
ed  at  positions  A  and  C,  which  are  identified 
by  an  angle  0,,  which  is  calculated  using  Eq. 
(5). 


The  first  step  in  the  simulation  is  the 
rotation  of  the  rigid  wheel  by  an  angle  u. 
Thus,  the  masses  move  to  locations  A  ’  and  C’ 
and  have  angles  0*  +  v.  In  the  second  step  the 
lumped  masses  deform  and  go  to  locations  A  ’  ’ 
and  C”.  Note  that  arbitrary  deformations 
where  chosen  to  illustrate  the  methodology 
used,  in  practice  the  masses  can  translate  in  all 
three  directions.  It  is  important  to  mention 
that  the  same  method  is  used  to  locate  the 
masses  in  each  of  the  three  planes. 

Even  though  the  lumped  masses  are  al¬ 
lowed  to  translate  in  all  directions,  it  is  as¬ 
sumed  that  the  displacements  are  very  small. 


such  that  a  translation  in  a  given  direction 
doesn’t  affect  those  in  other  directions.  This 
assumption  can  be  made  because  we  are  using 
a  very  large  number  of  masses.  This  is  a  valid 
assumption  for  this  particular  application  be¬ 
cause  this  model  is  to  be  run  on  relatively 
smooth  surfaces  at  slow  speeds,  which  guar¬ 
antees  that  we  will  not  be  getting  large  defor¬ 
mations  in  the  tire. 


Fig.  6:  Diagram  of  the  coordinate  systems  used  for  one 
plane.  Note  that  Y,y0  and  yt  are  positive  going  outside  of 
the  page 


1.4  Equations  of  motion 

For  illustrative  purposes  all  equations  present¬ 
ed  in  this  section  will  be  derived  for  Nm  is 
equal  to  1.  In  the  multibody  dynamics  com¬ 
munity  a  vector  of  generalized  coordinates  q 
represents  the  position  and  orientation  of  all 
bodies.  The  number  of  generalized  coordi¬ 
nates  for  this  model  is  9 Nm  +  6.  Thus,  for  the 
condition  that  Nm  is  equal  to  1,  vector  q  is 
defined  by  Eq.  (6),  where  i  is  equal  to  1.  Note 
that  if  Nm  is  greater  than  1  then  the  pattern 
will  continue;  left  sidewall  coordinates,  right 
sidewall  coordinates  and  finally  belt  and  tread 
coordinates.  It  is  important  to  mention  that  the 
first  index  (i)  represents  the  mass  number  and 
the  second  one  the  layer  (1-  left  sidewall,  r- 
right  sidewall,  bt-  belt  and  tread). 

q  = 

[x0  y0  z0  V  ip  ac  xui  yM  zui  (6) 

%i,r  yi,r  % i,r  ^i,bt  yi,bt  Z i,bt\ 
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Having  defined  the  generalized  coordinate 
vector  the  equations  of  motion  are  defined  by 
Eq.  (7), 

Mq  +  Cq  +  Kq  =  F  (7) 

where  the  vector  F  is  defined  by  the  applied 
external  forces,  C  and  K  are  the  damping  and 
stiffness  matrices,  and  M  is  the  generalized 
mass  matrix, 

'Mq  0  0  0  - 

0  Mu  0  0 

M  “  0  0  Mi  r  0  (8) 

_  0  0  0  Mi>bt_ 

1.1.1  Equations  of  motion  for  the  wheel 

The  mass  matrix  for  the  wheel  is  de¬ 


'  cos(0j  +  v)  0  sin(fy  +  v)' 

A2=  0  1  0  (13) 

sin  (6>i  +  u)  0  cos(0i  +  u)_ 

The  rotation  about  the  x-axis  from  the 
lumped  mass  reference  frame  to  the  global 
reference  frame  is  described  by  (14).  It  is  per¬ 
tinent  to  note  that  all  masses  have  the  same 
transformation. 

-10  0  - 
A3  =  0  cos  (ip)  sin(i p)  (14) 

.0  —sin  (ip)  cos(i p). 

The  rotation  about  the  z-axis  from  the 
lumped  mass  reference  frame  to  the  global 
reference  frame  is  also  the  same  for  all  masses 
and  it  is  described  by  equation  (15). 


fined  below, 

m0  0  0  0  0  0 

0  m0  0  0  0  0 

0  0  m0  0  0  0 

0  0  0  Iy  0  0 

0  0  0  0  Ix  0 

o  o  o  o  o  /z 

The  moments  of  inertia  ( Ix ,  lyand  lz ) 
are  calculated  in  Eq.  (10)  and  (11)  under  the 
assumption  that  the  wheel  is  a  solid  cylinder. 


4  =  4  = 


m0(3r02  +  b2) 


In  order  to  derive  the  equations  of  mo¬ 
tion  for  the  wheel  it  is  important  to  first  pre¬ 
sent  the  transformation  matrices  used  to  (i) 
change  between  reference  frames  and  (ii)  to 
orient  the  bodies  in  the  three-dimensional 
space.  The  transformation  matrix  from  the 
lumped  mass  reference  frame  to  the  global 
reference  frame  is  described  by  Eq.  (12),  not 
that  each  lumped  mass  has  a  different  trans¬ 
formation  matrix. 

"sin  (6>i  +  u)  0  —cos  (fy  +  v )' 

A1=  0  1  0  (12) 

Lcos(0j  +  u)  0  sin (0i  +  v )  . 

The  rotation  about  the  y-axis  for  each 
mass  has  the  following  transformation  from 
the  lumped  mass  reference  frame  to  the  global 
reference  frame, 


'  cos(ac)  sin(ac)  O' 

A4  =  —sin  (ac)  cos  (ac)  0  (15) 

0  0  1- 

The  forces  acting  on  the  wheel  due  to  the 
springs  and  dampers  are  found  in  the  Appen¬ 
dix  A,  Eq.(15) 

The  forces  due  to  the  springs  and 
dampers  acting  on  the  wheel  are  the  following 
in  the  global  reference  frame, 

4,o 

4.o 

4,o 

Observing  Eq.  (16)  it  is  evident  that  the 
forces  in  the  springs  and  dampers  connected 
to  the  wheel  are  multiplied  by  the  transfor¬ 
mation  matrices  to  bring  them  back  to  the 
global  reference  frame.  The  order  in  which 
these  matrices  are  multiplied  is  important.  In 
this  case,  the  tire  is  first  rotated  about  the  y- 
axis,  this  rotation  is  followed  by  a  rotation  of 
the  x-axis  and  finally  by  a  rotation  of  the  z- 
axis.  In  the  multibody  dynamics  community 
this  would  be  described  as  an  YXZ  rotation. 

The  equations  of  motion  for  the  wheel 
are  shown  in  Eq.  (17).  The  term  W  is  the  ap¬ 
plied  vertical  load  on  the  wheel  center,  while 
T  is  the  applied  driving  torque,  Mx  is  the  ap¬ 
plied  torque  about  the  x-axis,  and  Mz  is  the 
applied  torque  about  the  z-axis. 


Nm  [4/ 

^jA4a3a1  4,i 
i= 1  4  ,i 
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M0q0  — 


1  x,0 

Fy,  o 

Fz,o  ~  ™0g  -  W 
T 

Mx 

M7 


(17) 


The  following  equations  are  used  to 
obtain  the  forces  and  moments  at  the  wheel  in 
the  global  reference  frame, 

(18) 

(19) 

(20) 


Fx  —  Fx,o 

Fy  =  Fy,  0 

Fz  =  Fz0 
Nm 

=  X  U‘r° 

i= 1 


(21) 


Nr, 


Txp  =  ^  ~FU  sin(0j  +  v ) 

i=l 


yyt 


Nr, 


-Fu  cos  (6i  +  v)yyi 
~  F7,iT0  sign(sin(0i 

+  v)) 


(22) 


Tac  =  ^  -Fu  cos  (0t  +  v )  yyt 


i= l 


(23) 


The  equations  of  motion  for  the  left 
sidewall  are  defined  by  Eq.  (25).  The  first 
terms  are  due  to  the  internal  forces  in  the 
springs  and  dampers  attached  to  the  mass, 
while  the  second  set  of  terms  represent  the 
external  force.  The  ground  forces  will  be  ex¬ 
plained  in  detail  in  Section  2. 

The  equations  of  motion  for  the  right  sidewall 
are  presented  in  Eq.  (26) 


i,l 


Fx,i,l  +  F$i,i 

Fy,i,l 

Fz,i,l  T  ^3,uJ 

—mi  sg  cos (0j  +  u)  +  Fxcp ,  M 

—  Fycp,i,l  "E  Fyftd  i  i 

—misg  sin  (fy  +  v)  -  miRi:lv2  +  onXX 


(25) 


+  F2X  sin  (0j  +  v)  yyi 
-  F7Xr0sign( cos(0j 
+  u)) 

where  yy*  is  the  relative  location  of  the  mass 
in  the  lateral  direction  with  respect  to  the  belt 
and  tread  layer.  The  command  used  in 

Eq.  (22)  and  (23)  returns  a  -1  if  the  correspond¬ 
ing  element’s  angle  is  less  than  zero,  a  1  if  it 
is  greater  than  zero  and  a  0  if  it  equals  zero. 

1.1.2  Equations  of  motion  -  sidewalls 

The  mass  matrix  for  the  sidewalls  is  the 
for  both  sidewalls  since  it  is  assumed  that  an 
number  of  masses  is  used  in  each  sidewall. 

0 
0 


Fx,i,r  "E  F^i  r 
F  ■ 

1  y,i,r 

Fz,i,r  E  F3  i  r 


MixQi.r 


(26) 


mLSg  cos(6»i  -I-  u)  -I-  Fxcp X  r 
Fycp,i,r  —  FybdXir 

[~mi,sd  sin(0j  +  u)  —  mi  sRi  rv2  +  aUi 
Finally,  the  position  vectors  of  the 
lumped  masses  in  the  left  and  right  sidewall 
are  described  by  Eq.  (27)  and  (28). 


ri,i  = 


xQ 

To 

Lz0J 


R 


same 

equal 


I'i.r  ~ 


Mi, i  =  Mir  = 


mi>s  0 


m. 


771/ 


E  A4A3A2 

E  A4A3A2 

E  A4A3A2 


x0 

yo 

Lz0J 


0 

~zi,l 

Vi, i 

L  Xix  j 

Rs 
w 


(27) 


(24) 


0 

0  0 

Thus,  the  forces  acting  on  the  springs 
and  dampers  connected  to  the  left  sidewall 
(Eq.  (9-2),  (9-3)  and  (9-4))  and  right  sidewall 
(Eq.  (9-5),  (9-6)  and  (9-7))  are  presented  in 
Appendix  B. 


E  A4A3A2 


2 
0 

i,r 

Vl,r 

L  Xi,r  j 


(28) 


— z 


1.1.3  Equations  of  motion  -  belt  and  tread 

The  mass  matrix  for  the  belt  and  tread  layer  is 
defined  below, 
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Mi,bt  ~ 


mi,bt 

0 

0 


0 

^ i,bt 

0 


0 

0 

™i,bt 


(29) 


The  forces  acting  on  the  belt  and  tread  ele¬ 
ments  due  to  the  connected  springs  and  damp¬ 
ers  and  equation  of  motion  are  described  in 
Appendix  C,  Eq.  (9-8)  and  (9-9),  respectively. 

Using  the  same  formulation  as  the 
sidewall  masses,  the  location  of  the  lumped 
masses  in  the  belt  and  tread  layer  is  deter¬ 
mined  using  the  following  Eq.  (30) 


(30) 


-x0- 

Ru 

ri,bt  ~ 

yo 

+  A4A3A2 

0 

*0. 

.  0  . 

r — z 

i,bt 

+  a4a3a2 


yi,bt 
.  %i,bt  . 


2  Tire  -  soft  soil  interaction 

This  section  reviews  the  soil  model  used; 
thus,  it  goes  into  the  different  formulations 
used  to  calculate  the  shear  stresses  and  the 
normal  pressure.  Moreover,  it  goes  into  the 
details  of  how  the  interaction  between  the  tire 
and  the  soil  model  is  done.  A  general  sche¬ 
matic  of  the  different  variables  used  is  shown 


in  Fig.  7. 


Most  soils  have  elastic  and  plastic 
properties.  Thus,  the  ideal  behavior  of  an 
elasto-plastic  material  is  described  by  Fig.  10. 
However,  parameter  identification  for  the 


elastic  properties  of  soils  is  difficult.  As  such, 
the  theory  of  plasticity  will  be  used  to  model 
the  soil.  This  assumes  that  all  deformation  in 
the  soil  is  permanent,  so  the  rut  and  the  sink- 
age  are  assumed  to  be  the  same. 


✓ 

/ 

✓ 

/ 

_ 

Strain 


Elastic  Region 
Plastic  Region 


Fig.  8:  Ideal  behavior  of  an  elasto-plastic  material 


By  discretizing  the  tire  into  multiple 
masses  and  planes,  a  better  resolution  can  be 
obtained  for  the  forces  in  the  contact  patch. 
As  shown  in  Fig.  9,  the  width  of  the  tire  ( b )  is 
divided  into  three  sections  that  represent  the 
sidewall  and  tread  and  belt  layers.  The  formu¬ 
lation  for  the  width  of  each  section  (w,)  is  only 
dependent  on  the  width  of  the  tire, 


On  the  other  hand,  the  length  (si)  of  each  sec¬ 
tion  is  dependent  on  the  deformed  radius  and 
the  number  of  masses.  The  formulation  gives 
Eq.  (32). 

2n 

si  =  RiTT  (32) 


Fig.  9:  Discretization  of  the  contact  patch 

Before  the  formulation  for  the  ground 
forces  is  presented  it  is  pertinent  to  define  the 
longitudinal  slip  or  slip  ratio  (so).  The  slip  is 
important  when  it  comes  to  off-road  locomo¬ 
tion  because  without  any  slip,  motion  would 
not  be  possible.  The  slip  ratio  in  this  thesis  is 
calculated  using  the  SAE  standard,  which  is 
defined  by  the  equation  below, 
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s“=i-^  (33» 

where  Vx  is  the  longitudinal  velocity  of  the 

vehicle  and  co  is  the  radial  velocity  of  the 
wheel. 

2.1  Normal  Ground  Stress 

One  of  the  most  important  aspects  of  the 
interaction  of  the  tire  with  soft  soil  is  the 
computation  of  the  normal  pressure.  The  nor¬ 
mal  pressure  is  assumed  to  be  perpendicular 
to  the  surface  of  the  tire  contact  patch  and  it  is 
dependent  on  the  sinkage  of  the  tire;  a  dia¬ 
gram  of  this  contact  force  is  shown  in  Fig.  10. 
The  pressure  distribution  in  the  longitudinal 
direction  is  depicted  in  the  diagram  on  the  left 
of  Fig.  10.  At  the  entry  angle  the  pressure  in 
the  radial  direction  is  equal  to  zero  and  it  in¬ 
creases  until  the  location  of  peak  pressure, 
which  according  to  Reece  and  Wong  [2]  will 
shift  with  skidding  or  slipping;  the  pressure 
will  then  decrease  until  it  reaches  zero  at  the 
exit  angle.  On  the  other  hand,  the  pressure 
distribution  in  the  lateral  direction  is  symmet¬ 
ric  left  to  right,  as  shown  in  the  diagram  to  the 
right  of  Fig.  10.  However,  for  steering  ma¬ 
noeuvres  or  straight  line  driving  with  a  cam¬ 
ber  angle,  the  distribution  is  no  longer  sym¬ 
metric  about  the  centre.  The  shape  of  the  con¬ 
tact  patch  will  change  from  a  rectangle  to  a 
trapezoidal  shape. 

In  this  study  we  implement  an  approach 
developed  by  Grahn  [3,  4]  to  calculate  the 
normal  pressure  (pi),  which  is  a  method  de¬ 
rived  from  Bekker’s  formulation [5],  but  with 
additional  features;  it  includes  the  effects  of 
longitudinal  slip,  vehicle  velocity,  and  vertical 
penetration  velocity.  Grahn  found  that  pene¬ 
tration  velocity  has  a  large  influence  on  the 
pressure-sinkage  relationship  [6],  thus  motivat¬ 
ing  the  need  to  account  for  it  in  our  model. 


Fig.  10:  Left-  Longitudinal  pressure  distribution  in  the 
contact  patch.  Right  -  Lateral  pressure  distribution  in  the 
contact  patch  for  straight  line  driving.  The  red  arrow 
depicts  the  location  of  peak  pressure 

Using  this  methodology,  the  sinkage 
(Zj)  is  computed  for  each  section  i  in  the  dis¬ 
cretized  tire  model.  The  pressure  in  the  radial 
direction  on  each  mass  can  be  computed  using 
Eq.  (35).  Note  that  if  m  =  0  then  these  equa¬ 
tions  yield  Bekker’s  quasi-static  formulation 
[3].  Thus,  the  model  uses  Eq.  (34)  to  calculate 
the  pressure.  Furthermore,  it  is  important  to 
note  that  Eq.  (35)  is  only  valid  for  penetration 
velocities  larger  than  0.1  cm/s.  Thus,  when 
this  condition  is  not  met  the  model  uses  (34). 

Pi  =  (~s^  +  k(p) (34) 


Pi 


Ri  )l 


Once  the  pressure  is  obtained  at  each 
section,  the  ground  force  (<Jnp  is  calculated 
using  the  equation  shown  in  Eq.  (36), 

°n,i  =  PiWiSi  (36) 

As  it  was  previously  mentioned,  the  theory  of 
plasticity  is  used  to  model  the  soil.  Thus,  the 
ground  force  is  applied  to  the  tire  only  when 
the  vertical  velocity  of  the  wheel  is  either 
negative  (tire  sinks  into  the  ground)  or  equal 
to  zero  (tire  touches  the  surface). 

The  first  step  in  the  simulation  is  to 
load  the  tire.  Solving  the  equations  of  motion 
yields  the  position,  velocity  and  acceleration 
components  of  the  wheel  and  lumped  masses. 
Thus,  the  “find'  command  is  used  in 


arccos  1 


(35) 
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MATLAB  to  determine  the  first  mass  that  has 
sinkage,  which  in  turn  returns  the  entry  angle 
(de).  Similarly,  the  “finer  command  is  also 
used  to  determine  the  last  mass  that  has  sink- 
age,  yielding  the  exit  angle  (6r).  Consequent¬ 
ly,  higher  number  of  masses  ( Nm )  will  yield 
more  accurate  entry  and  exit  angle. 

The  second  step  in  the  simulation  is  to 
start  rolling  the  tire.  This  steps  starts  from  the 
initial  conditions  identified  in  the  dynamic 
settling  performed  in  step  one.  For  this  scenar¬ 
io  the  entry  angle  is  found  in  a  similar  man¬ 
ner.  However,  for  the  determination  of  the 
exit  angle,  the  “ max ”  command  is  used.  This 
command  will  find  the  mass  that  has  the  high¬ 
est  sinkage,  and  thus,  yield  the  exit  angle.  It  is 
important  to  note  that  an  algorithm  is  also 
used  to  determine  if  there  are  other  masses 
that  have  the  same  maximum  sinkage.  This 
guarantees  that  if  the  contact  patch  is  flat,  the 
correct  exit  angle  is  found. 

Having  computed  the  normal  pressure, 
the  total  drawbar  pull  (DP)  of  the  tire  is  the 
following, 

N_m 

DP  —  ^  '  ( Pxcp,i,l  @i,l 

i= 1 

—  On,i,l  COS  0M) 

N_m 

/  (  Pxcp,i,r  sift  @i,r 

t=i  (37) 

—  <Jn  i  r  COS  0;  r) 

N_m 

+  ^  {Pxcp  ,i,bt  sin  eUbt 
i= 1 

—  & n,i,bt  COS 

where  Fxcp  is  the  lateral  shear  force,  which  is 
referred  to  as  the  shear  trust  in  the 
terramechanics  community.  The  computation 
of  this  force  will  be  addressed  in  the  following 
section.  On  the  other  hand,  the  ground  force 
term  is  referred  to  as  the  compaction  re¬ 
sistance. 

In  contrast,  the  computation  of  forces 
in  the  vertical  direction  is  defined  by  Eq.  (38), 


—  ^ '  ( J'xcp.i.i  cos  On  +  i  i  sin  9iti) 
i= 1 

Nm 

3”  fxcp,i,r  COS  di,r 


i=l 

+  anir  sin  0i  r) 

Nm  ' 


y  (F  . 

/  \l  xcp,i, 


i,bt  COS  9i  bt 


3"  o? i,i,bt  sin 


where  W  is  the  normal  load  applied  to  the 
wheel. 


2.2  Shear  Stress 

The  shear  stress  allows  the  tire  to  have 
friction  with  the  ground  and  be  able  to  move. 
In  this  model  the  shear  stress  in  the  contact 
patch  is  computed  using  Janosi  and 
Hanamoto’s  approach  [7].  This  approach  uses 
the  tangential  velocity  or  interface  velocity  to 
calculate  the  shear  displacement  in  the  longi¬ 
tudinal  direction  (Jx  i)  for  each  mass  in  the 
contact  patch.  The  soil  shear  displacement 
will  be  zero  at  the  entry  angle  and  it  will  in¬ 
crease  until  it  reaches  the  exit  angle. 

ix,i  F-i  [(^i  9e)  (cos  Oe  cos  df)  ~9 

+  sd(cos6e  -  cosdi)]  t  J  ’ 
On  the  other  hand,  the  lateral  shear 
displacement  (jy  i)  is  calculated  using  the  lat¬ 
eral  slip  velocity.  The  lateral  soil  shear  dis¬ 
placement  behaves  in  the  same  way  as  in  the 
longitudinal  displacement  where  it  starts  at 
zero  and  keeps  increasing  until  it  reaches  the 
exit  angle. 

jy  i  Ri(  1  Sd)(9i  9e)  tan  ac  (40) 
Once  the  shear  displacements  are  calculated, 
the  longitudinal  and  lateral  shear  stress 
(tx  and  Ty)  can  be  evaluated  using  Eq.  (41) 
and  (42)  respectively.  This  approach  uses  the 
Mohr-Coulomb  soil  strength  failure  criterion. 

TXii  =  (c  +  ani  tan(0))  (l-e  k*  ^  (41) 
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T y,i  =  (c  +  anX  tan(0))  1  -  e 


hh 


(42) 


where  c  is  the  soil  cohesion,  0  is  the  internal 
angle  of  friction,  and  kx  and  ky  are  the  longitu¬ 
dinal  and  lateral  shear  deformation  modulus, 
respectively. 

Equations  (41)  and  (42)  are  valid  for 
either  pure  longitudinal  or  pure  lateral  simula¬ 
tions.  However,  for  combined  slip  scenarios 
the  shear  strength  envelope  of  the  soil  needs 
to  be  considered.  The  approach  followed  in 
the  thesis  is  consistent  with  that  developed  by 
Sandu  and  Chan  in  [8]  where  the  Mohr- 
Coulomb  soil  strength  failure  criterion  is  used 
to  define  the  maximum  shear  strength  of  the 
soil, 

T'max.i  =  (c  +  anX  tan (0))  (43) 

Thus,  the  maximum  shear  stress  envelope  is 
characterized  by  an  elliptical  constraint  de¬ 
fined  by  the  following  equation, 


max, i 
2 


+ 


Ly,i 


<  1 


(j-x,i)  T  {ji,max) 


e=e„ 


(44) 

(45) 


where  6a  is  the  angle  where  the  stresses  tran¬ 
sition  from  sticking  to  sliding.  Simplifying 
Eq.  (45)  yields  Eq.  (9-10),  which  can  then  be 
solved  numerically  to  find  9a,  see  Appendix 
D,  Eq.(9-10).  Consequently,  the  longitudinal 
and  lateral  shear  stresses  are  defined  by  the 
following  equations, 


LX,l 


f  ( 

ix,i\ 

Tmax,i  (  1  —  e 
<  ' 

kx 

9e>  9  >9a 

(46) 

k  ^ max,i 

Ty,i 

9a>  9  >  9r 

r  ( 

~hl\ 

^ max,i  1  1  —  ^ 

ky  I 

9e  >  9  >  9  a 

(47) 

v  T max,i 

9a>  9  >er 

In  principle,  Eq.  (46)  and  (47)  guarantee  that 
the  shear  strength  of  the  soil  never  exceeds  the 
failure  envelope,  however,  the  bulldozing 
force  is  not  included  in  this  formulation. 
Therefore,  there  might  be  times  were  the  shear 
strength  of  the  soil  could  be  exceeded  since 


the  bulldozing  is  added  on  top  of  the  shear 
force.  However,  this  simple  approach  is  a 
good  approximation,  since  the  bulldozing 
stress  is  relatively  small  compared  to  the  lat¬ 
eral  shear  stress. 

Using  the  same  approach  as  in  the  calculation 
of  the  normal  force,  the  longitudinal  and  lat¬ 
eral  shear  forces  acting  on  each  mass  are 
shown  in  Appendix  D  Eq.  (9-11)  and  (9-12), 
which  utilize  fix  and  ny  as  the  friction  coeffi¬ 
cients,  computed  using  Eq.  (48).  The  use  of 
these  coefficients  is  important  because  it  lim¬ 
its  the  available  traction  at  the  contact  patch  to 
1,  thus,  assuring  that  the  shear  forces  will  not 
exceed  the  shear  envelope.  Furthermore,  since 
the  shear  forces  are  bounded,  the  solver 
doesn’t  encounter  any  convergence  issues  due 
to  the  sliding  forces  exceeding  the  maximum 
shear  strength  of  the  soil. 

(Lmax.i  -  \  //io\ 

~  1  1  (48) 


2.3  Bulldozing  Force 


The  bulldozing  effect  is  created  when  a 
volume  of  the  soil  in  the  contact  patch  is  dis¬ 
placed  to  the  sidewall  of  the  tire  when  corner¬ 
ing.  As  such,  a  lateral  force  is  created  on  the 
sidewall  of  the  tire.  A  graphical  representation 
of  this  phenomenon  can  be  observed  in  Fig. 


Fig.  11:  Representation  of  the  bulldozing  force  [8] 

Few  tire  models  account  for  the  bull¬ 
dozing  component  of  the  lateral  force.  Those 
that  incorporate  it  applied  it  at  the  single  mass 
that  represents  the  contact  patch.  In  the  cur¬ 
rent  study,  due  to  the  discretization  of  the  tire, 
the  bulldozing  effect  can  be  applied  directly 
on  the  sidewall  elements,  thus  increasing  the 
realism  of  the  model.  However,  it  is  important 
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to  note  that  this  bulldozing  force  can  only  be 
applied  at  the  location  of  the  discrete  masses 
representing  the  sidewall. 

Thus,  a  bulldozing  force  will  be  ap¬ 
plied  at  each  mass,  which  is  also  a  more  real¬ 
istic  representation.  The  formulation  used  to 
determine  the  bulldozing  effect  is  based  on 
the  principle  of  passive  ground  resistance  de¬ 
veloped  by  Terzaghi,  which  is  presented  in  [9- 
11].  This  formulation  has  different  dimension¬ 
less  soil  resistance  coefficients  (Ny,  Nc,  Nq) 
that  according  to  Wong  [12]  are  dependent  on 
the  soil  angle  of  friction  0.  On  the  other  hand, 
q  is  the  surcharge  load  from  the  accumulated 
soil,  which  is  calculated  by  assuming  that  the 
soil  displaced  in  the  lateral  direction  on  each 
side  of  the  wheel  is  the  same  as  the  volume  of 
soil  displaced  by  the  area  shaded  in  grey  in 
Fig.  11.  Thus,  the  surcharge  load  q  is  discre¬ 
tized  by  the  method  used  in  Fig.  12  and  found 
using  the  equation  below. 


Where 


_  Ys^soil,i 
qi  mbZi 

(49) 

(n  0\ 
mb  «  tan  +  — J 

(50) 

Vsoil.i  ~  tan(tfc) 

(51) 

~  zisi 

(52) 

i  =  Ri(cos9e  —  cos  9i) 

(53) 

The  Terzaghi  soil  resistance  coefficients  [13] 
are  found  using  the  succeeding  equations, 

(54) 


N«  = 


2  cos 


(?4) 


iVc  =  c  010^-1)  (55) 

Ny  =  2{NC  +  1)  tan(0)  (56) 

Once  the  surcharge  is  determined  the  bulldoz¬ 
ing  force  is  obtained  using  Eq.  (57). 


Fybd,i  —  Ai(yszi  Ny  +  CZiNc 
+  q^iNy)  COS  0 


(57) 


It  is  important  to  note  that  the  bulldoz¬ 
ing  force  is  only  applied  to  the  outside  side- 
wall  when  cornering,  thus,  assuming  that  no 
soil  is  displaced  to  the  inside  of  the  corner. 
Furthermore,  it  is  also  important  to  mention 
that  it  is  estimated  that  Eq.  (57)  is  applicable 
to  slip  angles  under  45°. 


Top  View 


Fig.  12:  Volume  discretization  for  the  calculation  of  the 
bulldozing  force 

The  total  lateral  force  generated  by  the 
tire  is  the  summation  of  the  lateral  forces  at 
each  layer  plus  the  lateral  force  created  by  a 
camber  angle, 

Fy  —  Fy,i,l  T  Fy,i,r  T  Fy,i,bt  T  Fcamber 
Nm 

~  ^  '  (~Fycp,i,l  Fybd.i.l) 

i= 1 
N-m 

T  ^f^f—Fycp,i,r  —  Fybd,i,r )  (58) 

i=l 

Nm 

+  ^ ''f—Fycp.i.l ) 
i=l 

+  W  sin0 

2.4  Response  to  Repetitive  Loading 

The  multi-pass  is  an  important  effect  for 
soft  soil  interaction.  The  impact  of  a  tire  on 
the  soil  is  dependent  on  a  variety  of  factors. 
One  such  factor  is  the  type  of  pass.  A  towed 
wheel  will  not  alter  the  properties  of  the  ter¬ 
rain  in  the  same  manner  as  a  driven  wheel, 
which  induces  larger  changes  [2].  The  number 
of  passes  also  has  an  influence  on  the  terrain 
properties.  Moreover,  the  slip  ratio  also  af¬ 
fects  the  soil  properties  for  repetitive  loading 
scenarios. 

Even  though  other  researchers  such  as 
Harnisch  et  al.  [14]  have  used  Wong’s  ap¬ 
proach  with  satisfactory  results,  this  study  will 
rely  on  the  work  of  Holm  [15]  and  Senatore 
[16]  to  predict  the  response  to  the  multi-pass 
effect.  The  main  idea  of  this  approach  is  to 
estimate  new  parameters  for  the  soil  based  on 
the  pass  type,  slip,  and  number  of  passes. 
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Observing  the  results  of  Holm  in  [15]  it  is 
evident  that  the  drawbar  pull  increases  as  the 
number  of  passes  increases.  This  is  mainly 
due  to  the  fact  that  the  compaction  of  each 
successive  pass  reduces  the  rolling  resistance. 
Moreover,  it  is  also  evident  that  driven  wheels 
will  create  more  compaction  than  towed 
wheels,  and  thus,  have  higher  drawbar  pulls  as 
the  number  of  passes  increases.  Furthermore, 
it  also  evident  that  the  greatest  variation  in 
drawbar  pull  occurs  between  the  first  and  se¬ 
cond  pass.  Thus,  the  following  relationships 
were  developed  to  predict  the  change  in  soil 
conditions  due  to  repetitive  loading, 


C  —  (c  +  Csd )  +  C-^Sd 

(59) 

Ys  =  (Ys  +  Ysd )  +  FsiSd 

(60) 

kx  —  (Mx  T  ksd)  +  kxlsd 

(61) 

where  csd,  ysd,  and  ksd  are  parameters  that 
are  dependent  on  the  type  of  pass,  and  c^y^, 
and  kxl  are  functions  of  number  of  passes.  In 
the  case  of  a  towed  wheel,  csd,  ysd,  and  ksd 
are  equal  to  zero.  On  the  other  hand,  for  driv¬ 
en  wheels  their  value  is  dependent  on  the 
wheel  slip.  Instead,  for  the  parameters  that  are 
dependent  on  the  number  of  passes,  their  val¬ 
ue  is  equal  to  zero  for  the  first  pass  and  then  it 
is  contingent  to  the  number  of  passes. 


(235/85R16),  as  shown  in  Fig.  13.  According 
to  Michelin  of  North  America  [17]  this  is  an 
all  season  off-road  tire  for  a  light  truck.  It  has 
a  max  load  of  30421bs  (13,53  IN)  at  80psi 
(551.6kPa)  (load  index  120)  and  a  speed  rat¬ 
ing  of  R. 

0  Upgraded  tread  compound  for  tough  off-road  endurance  and  a  long 
tread  life. 

0  Three  circumferential  grooves  to  prevent  hydroplaning. 

0  Deep  tread  with  Michelin  biting  edges  for  exceptional  traction. 

0  Spiral  nylon  wrap  for  high  speed  endurance. 

0  Two  super  high  tensile  steel  belts 


Fig.  13:  Michelin  LTXA/T2  tire  used  for  the  model  valida¬ 
tion  [17] 

Due  to  the  fact  that  we  didn’t  have 
enough  test  data  for  this  tire,  the  parameters 
used  are  a  mixture  of  measured  and  approxi¬ 
mated  values.  All  damping  coefficients  were 
approximated  based  on  a  recommendation 
made  by  Mr.  Michael  Gipser  from  FTire; 
thus,  the  damping  coefficients  are  between  1- 
3%  of  the  stiffness  values  in  the  respective 
direction.  The  lateral  tire  stiffness  was  provid¬ 
ed  by  Dr.  Said  Taheri  as  an  approximation 
based  on  a  very  similar  tire.  The  sidewall,  belt 
and  tread  masses  are  estimations  based  on  the 
total  mass  of  the  tire.  The  remaining  values 
were  obtained  by  experimental  means,  such  as 
what  is  presented  in  Section  3.2.1. 


3  Case  Studies  and  Simulation  Results 

This  section  first  introduces  the  tire  and 
soil  parameters  used  for  all  the  simulations.  It 
then  presents  the  results  for  the  dynamic  load¬ 
ing  of  the  tire  on  rigid  ground.  The  experi¬ 
mental  results  are  compared  to  the  simulation 
results  in  an  effort  to  validate  the  tire  defor¬ 
mation.  The  succeeding  section  presents  a 
case  study  for  dynamic  settling  on  sandy 
loam.  Consequently,  case  studies  for  pure 
lateral,  pure  longitudinal  and  combined  slip 
scenarios  are  discussed.  Finally,  a  case  study 
for  repetitive  loading  is  analyzed. 

3.1  Tire  and  Soil  Parameters 

The  tire  chosen  for  the  validation  of  this 
model  is  Michelin  LTX  A/T  2  tire 


Table  2:  Tire  parameters  used  for  all  simulations. 
Note  that  all  this  parameters  need  to  be  divided  by  the 
number  of  masses  in  contact  and  mtotal  is  the  total  mass 
of  the  tire 


Parameter 

Value 

m  total 

19.83  kg 

mrim 

11.52  kg 

m  sidewall 

m,otai*(2/3)  kg 

mhelt  mtread 

m,0,ai*(l/3)  kg 

b 

0.2362  m 

Ru 

0.3937  m 

h 

0.0135  m 

ki 

215,172  N/m 

k2 

248,136  N/m 

k3 

10A5  N/m 

k4 

k2  N/m 

k5 

110,889  N/m 

k6 

5*10A5  N/m 

l<7 

110,889  N/m 
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ks 

10A6  N/m 

Cl 

ki*0.03  Ns/m 

C2 

k2*0.02  Ns/m 

C3 

k3*0.01  Ns/m 

c4 

k4*0.01  Ns/m 

C5 

ks*0.05  Ns/m 

C6 

k6*0.01  Ns/m 

C7 

k7*0.01  Ns/m 

The  sandy  loam  parameters  were 
found  from  a  source  written  by  Wong  [12], 
since  at  the  time  this  thesis  was  written  the 
extraction  of  the  soil  parameters  for  the  sandy 
loam  used  in  the  terramechanics  rig  was  still 
on-going.  It  is  important  to  note  that  the  index 
of  dynamic  pressure  could  not  be  found  either 
through  experimental  testing  or  literature; 
therefore,  the  Grahn  pressure-sinkage  equa¬ 
tion,  although  implemented  in  the  code,  can¬ 
not  be  used.  As  such,  all  the  simulations  pre¬ 
sented  in  this  paper  use  the  Bekker  formula¬ 
tion.  Grahn’ s  equation  can  be  employed  as 
soon  as  the  missing  parameters  are  experi¬ 
mentally  obtained. 

Table  3:Soil  Parameters  used  for  all  simulation 


Parameter 

Sandy  Loam 
(11%  Moisture 
Content) 

Cohesion  (c) 

4830  Pa 

Angle  of  friction 
(0) 

20° 

Index  of  static 
pressure  (n) 

0.9 

Coefficient  of  co¬ 
hesion  (kc) 

52,530  Pa 

Coefficient  of  fric¬ 
tion  (kp 

1,127,970  Pa 

Shear  stress  stiff¬ 
ness  (kx  &  ky) 

0.015  m 

Index  of  dynamic 
pressure  (m) 

_ 

Soil  Density  (y) 

1258  kg/m3 

3.2  Dynamic  Settling 

The  dynamic  settling  of  the  tire  is  an 
important  step  in  the  simulation  because  equi¬ 
librium  needs  to  be  attained  between  the  tire 
and  the  soil.  Moreover,  the  dynamic  settling 
of  the  tire  on  a  rigid  surface  was  used  to  vali¬ 
date  tire  deformation. 

3.2.1  Rigid  Surface 

The  first  step  in  the  validation  of  the 
model  was  to  compare  the  static  deflection  of 
the  tire  with  the  experimental  data.  The  data 
used  in  this  thesis  was  collected  by  Anake 
Umsrithong  from  Virginia  Tech  and  it  was 
collected  for  the  same  tire  used  in  the  soft  soil 
validation.  However,  it  is  important  to  note 
that  only  the  belt  deformation  will  be  validat¬ 
ed  since  no  data  is  available  for  the  sidewall 
deformation.  As  such,  for  tire  deflection  vali¬ 
dation  purposes,  the  model  was  run  using  only 
the  belt  and  tread  plane. 

The  setup  for  these  tests  was  simple.  It 
consisted  of  loading  and  unloading  the  tire  at 
different  loads  and  then  recording  the  deflec¬ 
tion  at  32  evenly  spaced  locations  around  the 
circumference  of  the  tire.  A  picture  of  the  test 
setup  can  be  observed  in  Fig.  14. 


Fig.  14:  Quarter  car  rig  testing  setup  used  to  validate  static 
deflection 

During  the  loading  of  the  tire,  pictures 
were  taken  perpendicular  to  the  tire  surface 
such  that  image  processing  could  be  used  to 
identify  the  deflections  in  the  tire.  Thus,  two 
different  loads  were  used  to  validate  the  mod¬ 
el.  The  first  one  is  an  intermediate  load  of 
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5,128  N  and  the  second  is  a  high  load  of 
1 1,022  N. 

The  comparison  between  the  experi¬ 
mental  and  simulation  results  for  the  interme¬ 
diate  load  can  be  observed  in  Fig.  15.  Each 
circle  represents  the  location  of  the  masses  in 
the  simulation  and  the  points  marked  in  the 
tire.  The  deformation  predicted  by  the  simula¬ 
tion  closely  predicts  the  deformation  observed 
in  the  experimental  testing. 


Fig.  15:  Static  deflection  comparison  for  a  rigid  surface 
under  an  intermediate  load  (5128N) 


Fig.  16  shows  the  deformation  of  the  lumped 
masses  in  contact  with  the  ground  for  both  the 
radial  and  circumferential  directions  when  a 
medium  load  (5,128  N)  is  applied.  By  analyz¬ 
ing  the  data  it  is  evident  that  there  is  only  a 
small  deformation  in  the  circumferential  di¬ 
rection,  while  the  largest  radial  deflection  is 
around  2.5  cm. 


Time(s) 


Time  (3) 


Time(s) 


0.01 


Time  (3) 


Fig.  16:  Lumped  mass  deformations  in  the  circumferen¬ 
tial  (x)  and  radial  (z)  direction  for  masses  in  contact 
(8,9,10)  at  a  medium  load  (5128N) 

The  results  for  the  high  load  also  show  good 
agreement  with  the  experimental  data  as  illus¬ 
trated  in  Fig.  17.  Moreover,  by  observing  Fig. 
17  and  Fig.  18  it  is  evident  that  the  defor¬ 
mation  experienced  by  the  tire  is  much  larger 
for  this  higher  load.  Thus,  the  close  resem¬ 
blance  between  the  collected  data  and  the 
simulation  suggests  that  the  model  is  able  to 
accurately  predict  tire  deformation  for  static 
conditions. 
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Fig.  17:  Static  deflection  comparison  for  a  rigid  surface 
under  high  load  (11,022  N) 
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Fig.  18:  Simulation  comparison  of  tire  deformation  for 
different  loads  under  a  rigid  surface 

3.2.2  Soft  Soil 

This  section  will  review  the  dynamic  loading 
simulations  that  were  done  on  sandy  loam. 
The  parameters  for  the  sandy  loam  soil  have 
already  been  presented  in  Table  3.  Thus,  the 
results  presented  in  this  section  and  the  suc¬ 
ceeding  sections  are  more  qualitative  in  na¬ 
ture.  When  available,  the  simulation  will  be 
compared  to  data  collected  for  the  same  tire  in 
the  AVDL  Terramechanics  rig,  which  will  be 
thoroughly  explained  in  Section  4.1.  Even 
though  the  sandy  loam  used  in  the  simulation 
most  closely  matches  the  silky  sand  in  the 
terramechanics  rig,  the  soils  are  not  identical. 
Therefore  these  comparisons  should  not  be 
taken  as  a  direct  validation,  but  more  of  as  an 
indication  that  the  tire  model  is  properly  pre¬ 
dicting  tire  behavior  for  this  type  of  soil. 
Moreover,  as  mentioned  previously,  due  to 
lack  of  soil  parameters,  the  pressure-sinkage 


equation  formulated  by  Grahn  couldn’t  be 
used,  as  such,  all  the  simulations  in  this  inves¬ 
tigation  use  Bekker’s  pressure-sinkage  equa¬ 
tion. 

The  dynamic  loading  of  the  tire  is  important 
because  it  is  key  to  the  analysis  while  equilib¬ 
rium  is  reached  between  the  tire  and  the 
ground.  Thus,  it  is  always  the  first  step  of  the 
simulation.  Fig.  19  shows  a  plot  of  the  vertical 
displacement  and  velocity  of  the  wheel  center 
as  the  tire  is  dropped  into  the  soft  ground.  As 
expected,  it  can  be  observed  that  the  wheel 
center  has  a  larger  vertical  displacement  as  the 
normal  load  is  increased.  Similarly,  the  verti¬ 
cal  velocity  follows  the  same  behavior. 


Time  (s)  Time  (s) 


Fig.  19:  Left-Vertical  displacement;  Right-  Vertical 
velocity  for  the  rigid  wheel  for  the  dynamic  settling  sce¬ 
nario 

Fig.  20  shows  the  normal  load  fluctuation  on 
the  wheel  center  as  the  tire  is  dropped  in  the 
soft  soil.  As  expected,  the  tire  bounces  slight¬ 
ly  until  it  reaches  equilibrium.  However,  by 
careful  inspection  is  it  evident  that  higher 
loads  cause  higher  fluctuations.  Moreover, 
due  to  the  increased  vertical  velocity  it  is  also 
observed  that  at  higher  loads  the  soil  experi¬ 
ences  higher  pressure.  It  is  important  to  note 
that  different  fluctuations  were  observed  for 
different  damping  coefficients,  however,  since 
this  thesis  is  not  interested  in  the  transient 
behavior  of  the  tire  this  is  not  of  great  conse¬ 
quence. 
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Time(s) 

Fig.  20:  Forces  at  the  wheel  centre  for  the  dynamic 
loading  scenario 

Besides  the  displacements  and  forces  at  the 
wheel  center  the  model  is  able  to  predict  rut 
height  for  both  the  sidewall  and  tread  and  belt 
layers.  It  is  important  to  note  that  since  a  plas¬ 
tic  model  was  adopted  for  the  soil,  the  rut 
height  will  be  equal  to  the  sinkage  of  the  tire. 
Observing  the  sinkage  of  both  the  sidewalls 
and  tread  and  belt  layers  in  Fig.  21  it  is  evi¬ 
dent  that  the  belt  and  tread  layer  has  a  larger 
sinkage  due  to  its  larger  radius.  In  Fig.  21  the 
simulation  results  were  compared  to  experi¬ 
mental  data  gathered  in  the  terramechanics  rig 
using  the  silky  sand  discussed  earlier. 


Fig.  21:  Sinkage  as  a  function  of  normal  load  for  the 
dynamic  loading  case.  Note  that  the  simulation  results  are 
obtained  using  a  Sandy  Loam  (11%  moisture  content), 
while  the  experimental  data  is  obtained  using  a  Silky  Sand 
(~2%  moisture  content ) 

Fig.  22  shows  the  pressure  distribution  in  the 
tire  for  a  normal  load  of  6,000  N.  The  tire  is 
divided  into  three  layers;  the  middle  layer 
represents  the  tread  and  belt  layer,  and  the  left 
and  right  layers  the  sidewalls.  Each  rectangle 
represents  a  lumped  mass,  thus,  it  is  assumed 


that  the  pressure  is  the  same  in  any  location  of 
that  specific  rectangle.  Moreover,  the  color 
bar  represents  the  normal  pressure  exerted  on 
each  lumped  mass  in  Pascal’s.  By  inspection, 
it  is  evident  that  the  pressure  distribution  is 
symmetrical  left-to-right,  however,  we  can  see 
that  the  pressure  is  higher  for  the  tread  and 
belt  layer,  which  makes  sense  since  it  has  a 
larger  sinkage. 


Fig.  22:  Pressure  distribution  for  the  dynamic  loading 
case  using  a  normal  load  of  6,000  N,  note  that  the  color 
bar  represents  the  normal  pressure  in  Pascal ’s. 

Finally,  in  Fig.  23  a  diagram  of  the  deformed 
tire  is  showed  for  the  6,000  N  case.  It  can  be 
observed  that  the  tire  experiences  a  considera¬ 
ble  amount  of  deformation  in  the  contact 
patch,  which  obviously  increases  with  an  in¬ 
crease  in  normal  load. 
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Fig.  23:  Deformed  tire  for  the  dynamic  loading  case 

using  a  normal  load  of  6, OOON 


3.3  Steady  State  -  Pure  Longitudinal  Slip 

The  pure  longitudinal  slip  case  study  was  run 
using  the  sandy  loam.  For  all  simulations,  a 
longitudinal  velocity  of  2  m/s  was  inputted, 
thus,  the  results  presented  in  this  section  show 
the  steady-state  behavior  of  the  tire.  Further¬ 
more,  72  lumped  masses  per  layer  of  the  tire 
model  were  used  for  all  simulations. 

Tire  sinkage  is  important  when  predicting  the 
performance  of  off-road  vehicles.  Fig.  24 
shows  the  simulation  results  for  sidewall 
sinkage  at  a  normal  load  of  2,000  N  for  differ¬ 
ent  slips.  At  first  inspection  it  is  evident  that 
there  is  a  lot  of  variation  in  the  data,  however, 
this  is  due  to  the  fact  that  masses  lose  and 
gain  contact  with  the  ground  constantly.  Ob¬ 
viously,  if  the  number  of  masses  is  increased 
it  is  expected  that  these  fluctuations  will  de¬ 
crease. 
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Fig.  24:  Sidewall  sinkage  for  different  longitudinal  slips 
at  a  normal  load  of 2,000  N  and  zero  slip  angle 


In  order  to  process  the  data,  mean  sinkage 
values  were  calculated  for  the  different  nor¬ 
mal  loads  and  slip  ratios,  which  are  presented 
in  Fig.  25.  In  the  same  way  as  in  the  dynamic 
loading  of  the  tire,  the  belt  and  tread  layer  has 
a  higher  sinkage  due  to  its  larger  radius.  It  is 
also  evident  that  sinkage  has  a  slight  increase 
for  slip  ratios  of  l%-60%.  However,  for  high¬ 
er  slip  ratios  higher  sinkage  is  observed;  this 
is  referred  as  slip  sinkage.  Essentially,  the 
slipping  of  the  tire  digs  in  the  ground  creating 
higher  sinkage  values. 


— ©— Side-2000N 
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— O — Side^IOOON 
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Side-6  00  ON 
— * —  B+T-6000N 


Longitudinal  Slip(%) 


Fig.  25:  Sinkage  as  a  function  of  slip  for  different  normal 
loads  at  zero  slip  angle 


The  drawbar  pull  is  probably  the  most  im¬ 
portant  measure  when  it  comes  to  off-road 
locomotion;  since  it  predicts  vehicle  mobility. 
Fig.  26  shows  the  drawbar  pull  as  a  function 
of  time  for  different  slip  ratios  at  a  normal 
load  of  2,000  N.  Similarly  to  Fig.  24,  fluctua- 
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tions  are  observed  in  the  data  due  to  the  loss 
and  gain  of  contact  of  masses  in  the  contact 
patch. 


Fig.  26:  Drawbar  pull  for  different  slip  ratios  at  a  nor¬ 
mal  load  of 2,000  N 

As  it  was  previously  mentioned,  an  in¬ 
crease  in  lumped  masses  will  decrease  fluc¬ 
tuations,  and  this  can  be  observed  in  Fig.  27. 
This  simulation  was  run  using  200  masses, 
rather,  than  the  standard  72  masses  used  for 
all  other  simulations.  The  fluctuations  are 
significantly  smaller  for  this  simulation;  how¬ 
ever,  the  mean  drawbar  pull  for  the  data  set  is 
quite  close  to  the  simulation  with  72  masses. 


Fig.  27:  Drawbar  pull  as  a  function  of  time  for  200  mass¬ 
es,  a  normal  load  of 2,000  N,  and  a  slip  ratio  of  30%. 

In  order  to  process  the  data  presented 
in  Fig.  26,  mean  averages  were  calculated  for 
each  normal  load  and  slip  ratio.  The  result  is 
Fig.  28.  As  expected,  an  increase  in  normal 
load  yields  higher  drawbar  pulls.  This  is  due 
to  the  fact  that  high  loads  create  higher  sink- 
age,  thus,  increasing  the  ground  pressure  and 


the  longitudinal  shear  stress.  Moreover,  it  can 
be  seen  in  Fig.  28  that  the  drawbar  pull  in¬ 
creases  with  longitudinal  slip  until  it  reaches 
40%,  where  it  starts  to  decrease.  This  transi¬ 
tion  point  defines  the  location  where  rolling 
resistance  stops  being  the  limiting  factor,  and 
the  wheel  slip  becomes  the  limiting  factor. 
Thus,  the  location  of  this  transition  point  is 
heavily  based  on  soil  properties. 


Fig.  28:  Drawbar  pull  as  a  function  of  longitudinal  slip 
for  different  normal  loads 


Comparing  the  simulation  results  of  Fig.  28 
with  experimental  data  collected  at  AVDL  in 
Fig.  29  for  a  rigid  wheel  on  GRC1  lunar  soil 
simulant  it  is  evident  that  the  data  follows  a 
similar  trend. 


32x11  Rigid  Wheel  DP  vs.  Sd 
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Fig.  29:  Drawbar  pull  vs.  slip  ratio  for  experimental 
results  on  a  rigid  wheel  collected  at  AVDL  using  GRC1 
lunar  soil  simulant  [18] 

The  driving  torque  is  another  im¬ 
portant  measure  when  it  comes  to  off-road 
mobility.  Knowing  the  required  torque  to 
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traverse  a  given  terrain  is  important  when 
determining  the  required  power  needed  for  a 
vehicle.  Fig.  30  shows  the  driving  torque  re¬ 
quired  to  traverse  the  terrain  at  various  normal 
loads  and  slip  ratios.  As  expected,  higher 
normal  load  yields  higher  torque. 


Fig.  30:  Driving  torque  as  a  function  of  longitudinal  slip 
for  different  normal  loads 

The  camber  angle  has  an  important  in¬ 
fluence  in  the  prediction  of  soft  soil  locomo¬ 
tion.  Most  importantly,  the  camber  has  a  large 
influence  on  the  pressure  distribution  in  the 
contact  patch.  Fig.  3 1  shows  the  pressure  dis¬ 
tribution  for  a  tire  travelling  at  zero  slip  angle 
and  zero  camber.  Similarly  to  Fig.  22,  it  can 
be  observed  that  the  normal  pressure  is  sym¬ 
metric  left  to  right.  On  the  other  hand,  Fig.  32 
shows  the  pressure  distribution  for  a  tire  trav¬ 
elling  at  2°  camber  and  zero  slip  angle.  In  this 
case,  higher  pressures  are  experienced  by  the 
right  sidewall. 


Fig.  31:  Pressure  distribution  for  a  normal  load  of  2,000 
N  and  a  longitudinal  slip  of  10%  at  zero  slip  angle,  the 
color  bar  represents  the  normal  pressure  in  Pascal ’s 


Fig.  32:  Pressure  distribution  for  a  normal  load  of  2,000 
N,  a  longitudinal  slip  of  10%,  zero  slip  angle  and  a  cam¬ 
ber  angle  of  2°,  the  color  bar  represents  the  normal  pres¬ 
sure  in  Pascal ’s 

Furthermore,  the  camber  angle  creates 
additional  lateral  forces  that  contribute  to  the 
overturning  moment.  Fig.  33  shows  both  the 
driving  torque  and  the  overturning  moment 
for  a  tire  traveling  with  a  camber  angle.  The 
driving  torque  is  not  affected  by  the  camber 
angle;  however,  a  small  overturning  moment 
can  be  observed.  As  it  was  previously  men¬ 
tioned,  the  fluctuations  in  the  data  are  due  to 
the  change  of  masses  in  contact. 
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Fig.  33:  Left-  Driving  torque  and  Right-  overturning 
moment  for  a  normal  load  of  2,000  N,  a  longitudinal  slip 
of  10%  and  a  camber  angle  of  2  ° 


3.4  Steady  State  -  Pure  Lateral  Slip 

The  next  case  that  will  be  investigated  is 
the  pure  lateral  slip.  Fig.  34  shows  the  lateral 
force  as  a  function  of  slip  angle  for  different 
normal  loads.  As  expected,  the  lateral  force 
increases  with  both  slip  angle  and  normal 
load.  As  it  was  shown  in  the  previous  section, 
higher  normal  load  yields  higher  sinkage, 
thus,  yielding  higher  lateral  shear  force. 
Likewise,  the  higher  sinkage  and  slip  angle 
will  also  increase  the  volume  of  soil  dis¬ 
placed,  consequently,  increasing  the  bulldoz¬ 
ing  force.  However,  similarly  to  the  drawbar 
pull  plot  in  Fig.  28  there  is  a  transition  point 
where  the  forces  reach  the  soils  shear  strength, 
thus,  they  level  out. 


Fig.  34:  Lateral  force  as  a  function  of  slip  angle  for 
different  normal  loads  at  zero  longitudinal  slip 

Fig.  35  presents  the  overturning  and 
self-aligning  moments  as  a  function  of  slip 
angle  for  different  normal  loads.  It  is  evident 
that  both  moments  follow  similar  patterns; 
both  increase  with  slip  angle  and  normal  load. 
Contrary  to  the  lateral  force,  the  moments 
start  leveling  out  at  around  25°  slip  angle. 


Slip  Angle  (deg) 


Fig.  35:  Left-  Overturning  moment;  Right-  Self-aligning 
moment  as  a  function  of  slip  angle  for  different  normal 
loads  at  zero  longitudinal  slip 


Finally,  Fig.  36  shows  the  pressure 
distribution  for  a  steered  tire.  Due  to  the  fact 
that  camber  is  kept  constant  at  0°,  the  change 
in  pressure  from  left  to  right  is  not  that  signif¬ 
icant.  However,  some  significant  deformation 
due  to  the  lateral  shear  and  bulldozing  forces 
is  encountered  in  the  lateral  direction  of  the 
contact  patch. 
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Fig.  36:  Pressure  distribution  in  the  contact  patch  for  a 
slip  angle  of  15°,  a  slip  ratio  ofO%  and  a  normal  load  of 
2,000  N.  Note  the  deformation  in  the  lateral  direction  of 
the  contact  patch 

3.5  Combined  Slip 

Pure  lateral  and  pure  longitudinal  case 
studies  are  quite  common  in  tire  modelling, 
however,  vehicles  are  usually  travel  under 
combined  slip  scenario;  therefore,  it  is  im¬ 
portant  to  account  for  combined  slip  in  the 
model. 

Fig.  37  shows  the  variation  in  drawbar 
pull  for  different  slip  angles  at  a  normal  load 
of  2,000  N.  The  figure  shows  that  an  increase 
in  slip  angle  reduces  the  drawbar  pull,  which 
is  consistent  with  the  observations  made  by 
various  researchers.  This  reduction  is  attribut¬ 
ed  to  the  fact  that  an  increase  in  slip  angle  will 
increase  the  lateral  shear  stress,  thus,  reducing 
the  longitudinal  shear  stress;  since  the  com¬ 
bined  forces  cannot  exceed  the  shear  strength 
of  the  soil.  However,  it  is  important  to  note 
that  at  negative  slip  ratios  the  drawbar  pull  is 
higher  than  at  their  respective  positive  slips. 
This  phenomena  is  attributed  to  the  fact  that 
braking  will  create  higher  sinkage,  thus,  in¬ 
creasing  the  longitudinal  shear  stress. 


Fig.  37:  Drawbar  pull  as  a  function  of  longitudinal  slip 
for  different  slip  angles  at  a  normal  load  of 2,000  N 


Fig.  38  shows  the  lateral  force  for  dif¬ 
ferent  slip  ratios  at  a  normal  load  of  2,000  N. 
As  anticipated,  the  lateral  force  decreases  with 
an  increase  in  slip  ratio.  Given  the  relationship 
between  the  longitudinal  and  the  lateral  shear 
stress,  an  increase  in  slip  ratio  will  increase 
the  longitudinal  shear  stress,  thus  reducing  the 
lateral  shear  stress  and  vice  versa. 
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Fig.  38:  Lateral  force  as  a  function  of  slip  angle  for 
different  longitudinal  slips  at  a  normal  load  of 2,000  N 

Fig.  39  shows  the  stress  distribution  at 
the  contact  patch  for  the  combined  slip  sce¬ 
nario.  As  it  can  be  observed,  the  longitudinal 
shear  stress  increases  very  quickly  from  the 
entry  angle.  However,  at  an  angle  9a,  the  tire 
starts  to  slide,  which  creates  a  drop  in  shear 
stress.  Similarly,  the  lateral  shear  stress  also 
increases  until  the  angle  9a,  where  it  also 
drops.  This  decrease  in  shear  stress  at  the  an¬ 
gle  9a  guarantees  that  the  Mohr-Coulomb 
failure  criteria  is  upheld 
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Fig.  39:  Shear  stress  distribution  for  the  belt  and  tread 
layer  at  a  normal  load  of  4, OOON,  a  slip  angle  of  10°  and 
40%  slip  ratio 


4  Experimental  Test  Work 

An  exceptional  amount  of  experimental 
validation  work  is  required  to  ensure  accurate 
and  robust  performance  of  the  newly  devel¬ 
oped  tire  model.  Therefore  the  following  sec¬ 
tions  will  highlight  the  indoor  test  rig  where 
the  main  source  of  data  collection  is  obtained, 
the  tire  test  subject,  additional  instrumenta¬ 
tion,  test  soil  and  laboratory  soil  parameter 
extraction. 


4.1  Indoor  Test  Rig 


3.6  Repetitive  Loading 

The  following  section  will  present  the  results 
for  repetitive  loading.  Again,  due  to  the  una¬ 
vailability  of  experimental  data  for  repetitive 
loading,  the  results  presented  are  more  quali¬ 
tative  in  nature.  Fig.  40  shows  the  drawbar 
pull  for  a  normal  load  of  2,000  N  and  a  longi¬ 
tudinal  slip  of  20%  for  three  consecutive  driv¬ 
en  tires.  The  simulation  predicts  that  drawbar 
pull  increases  with  each  pass,  which  is  due  to 
the  fact  that  less  sinkage  is  experienced,  thus, 
less  rolling  resistance.  This  also  agrees  with 
the  experimental  results  collected  by  Holm  in 
[15]. 


The  main  test  platform  for  investigating 
the  tire-soil  interaction  for  model  validation  is 
the  Terramechanics  Rig  in  the  Advanced  Ve¬ 
hicle  Dynamics  Lab  at  Virginia  Tech.  This 
single-wheeled  test  rig  essentially  simulates  a 
rigidized-suspension,  quarter-car  model  made 
to  extensively  investigate  tire  dynamics  when 
applied  to  a  wide  variety  of  terrains,  as  shown 
in  Fig.  41.  The  rig  is  driven  by  two  motors  - 
one  provides  the  longitudinal  velocity  of  the 
test  wheel  by  actuating  the  belt  system  which 
moves  the  carriage  (pictured  in  yellow  in  Fig. 
41)  along  the  test  terrain,  and  the  second  mo¬ 
tor  controls  the  angular  velocity  of  the  test 
wheel;  thus,  allowing  the  user  to  predefine  the 
different  longitudinal  slip  test  scenarios.  The 
rig  takes  measurements  via  the  Kistler  P650 
wheel  hub  sensor,  which  measures  all  three 
forces  and  moments  attributed  to  the  tire  -soil 
interaction,  wheel  rotation,  and  angular  veloc¬ 
ity,  also  a  string  potentiometer  is  used  to 
measure  the  wheel  vertical  travel.  More  in¬ 
formation  of  the  development  and  capabilities 
of  the  Terramechanics  rig  can  be  found  in  [19, 
20], 


Fig.  40.  Drawbar  pull  as  a  function  of  time  for  repetitive 
loading  using  a  normal  load  of  2,000  N,  zero  slip  angle 
and  a  longitudinal  slip  of  20% 
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Fig.  41:  Terramechanics  rig 


4.2  Design  of  Experiments 

A  design  of  experiments  (DOE)  test  ma¬ 
trix  has  been  developed  to  efficiently  accom¬ 
modate  the  needs  of  the  model  validation  and 
stay  within  the  capabilities  of  the  Terra- 
mechanics  rig,  as  shown  in  Table  4.  The  DOE 
will  be  performed  in  a  two  stage,  fraction  fac¬ 
torial  design  to  insure  statistically  relevant  and 
repeatable  results.  The  DOE  was  established 
with  the  following  main  interests  in  mind: 
high  slip  scenarios,  similar  parameters  of  an 
actual  full-scale  test  vehicle,  and  practicality 
of  soil  conditioning.  When  developing  the 
DOE,  the  military  advisor  of  the  project  ex¬ 
pressed  the  importance  and  interest  to  the 
military  to  test  at  large  slip  values.  The  actual 
test  tire,  the  Michelin  LT  235/85R-16  LTX 
A/T  2  E,  as  mentioned  in  Section  3.1,  was 
chosen  to  correlate  test  results  measured  in  the 
Terramechanics  test  rig  and  further  validate 
the  tire  model  with  the  results  of  an  instru¬ 
mented  test  vehicle  at  the  University  of  Preto¬ 
ria.  For  this  reason,  the  normal  load  and  infla¬ 
tion  pressures  were  established  to  mimic  that 
of  the  instrumented  Land  Rover  Defender 
110.  The  soil  parameters  were  established  to 
keep  the  soil  preparation  process,  elaborated 
in  section  3.2.9,  repeatable  and  relevant,  yet, 
less  strenuous  on  the  investigator. 


Table  4:  Design  of  experiments  test  matrix 


Parameter 

No.  of 
Levels 

Range 

Slip 

Percentage 

8 

20  -  100% 

Normal  load 

2 

5000  N  & 
6000  N 

Inflation 

Pressure 

2 

29  psi  & 
35%  psi 

Moisture 

Content 

2 

2-3%  &  ~9% 

Compaction 
Resistance/ 
Cone  Index 

2 

10.1  Kpa/mm 
& 

-11.1  kPa/mm 

Note,  aside  from  the  slip  percentage,  the  other 
parameters  have  an  arbitrary  second  level  for 
further  model  validation. 

4.3  Preliminary  Test  Results 

Following  the  design  of  experiments,  sets  of 
data  have  been  collected  to  reflect  the  lower 
level  parameters  of  the  test  matrix,  as  shown 
in  Fig.  42.  These  results  are  presented  as  nor¬ 
malized  drawbar  pull  versus  slip,  for  the  low 
level  parameters  of  the  DOE.  By  inspection 
of  Fig.  42,  one  can  see  that  the  drawbar  pull 
follows  the  expected  theoretical  relationship 
with  a  change  in  slip. 


slip  (%) 

Fig.  42:  Drawbar  Pull  Coefficient  vs.  slip  at  low  DOE  pa¬ 
rameters 

From  this  test  work,  it  became  evident  that 
some  equipment  of  the  rig  needed  revision 
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and  additional  instrumentation,  primarily  the 
normal  load  and  slip  control.  Throughout  the 
testing,  the  normal  load  applied  to  the  wheels 
varied  too  frequently  to  be  able  to  directly 
validate  the  model  performance  with  the  col¬ 
lected  data.  Therefore  in  the  following  subsec¬ 
tions,  new  instrumentation  to  control  wheel 
normal  and  improve  slip  calculation  will  be 
discussed. 

4.4  New  Instrumentation 

In  order  to  accurately  set  and  measure 
experimental  parameters  in  the  Terramechan- 
ics  rig,  it  is  vital  to  include  more  instrumenta¬ 
tion.  The  new  instrumentation  and  their  pur¬ 
pose  will  be  discussed  in  the  following  sub¬ 
sections. 

3.2.3  Measure  Deformed  Rolling  Radius 

As  mentioned  in  the  previous  section, 
the  Terramechanics  rig  can  run  a  predefined 
slip,  assuming  the  user  knows  the  rolling  de¬ 
formed  radius.  While  this  is  normally  deter¬ 
mined  through  repetitive  measures  of  the 
wheel’s  circumference  after  rolling  at  an  as¬ 
sumed  0%  slip  condition,  error  is  expected 
due  to  the  nature  of  needing  a  deformed  roll¬ 
ing  radius  to  calculate  slip.  Therefore,  a 
method  to  accurately  measure  the  rolling  de¬ 
formed  radius  is  taken  upon  the  rig,  which 
also  gives  an  opportunity  to  measure  the 
wheel  contact  patch  deformation  and  max 
sinkage  in  soil. 

The  approach  taken  to  measure  the  roll¬ 
ing  deformed  radius  of  a  tire  is  to  place  multi¬ 
ple  optical  distance  measuring  sensors  inside 
of  a  tire  to  measure  the  deflection  of  the  tire 
inner  liner.  Multiple  strategies  and  sensors 
have  been  used  in  the  past  to  instrument  tires 
to  collect  deformation  data.  Sensor  technology 
used  to  instrument  tires  range  among  several, 
such  as  optical  [21-25],  ultrasonic  sound  [26], 
passive  surface  acoustic  wave  (SAW)  [27], 
magnetic  and  piezoelectric[28-31]. 

The  essential  components  that  make 
the  Wireless  Internal  Tire  Sensors  (WITS) 
possible  are  highlighted  as  the  data  collection, 
wireless  data  transmission,  and  power  supply. 


Starting  with  the  data  collection,  multiple 
SHARP  infrared  emitting  diodes  (IREDs)  are 
placed  along  different  angular  positions  of  the 
rim  to  collect  as  much  tire  deflection  data  as 
possible,  as  shown  in  Fig.  43.  Given  the 
small  test  path  of  the  Terramechanics  rig 
(roughly  4  meters),  the  test  wheel  commonly 
makes  only  two  full  rotations  between  testing 
so  the  more  sensors,  the  more  data.  Also,  due 
to  the  low  speed  testing  conditions  of  the  tire, 
normally  5-20  cm/s,  tire  imbalance  is  consid¬ 
ered  negligible;  despite  the  assumption,  ef¬ 
forts  to  balance  the  sensor  placement  is  al¬ 
ways  taken  into  consideration. 


Fig.  43:  Sample  diagram  of  sensor  arrangement  inside  a 
wheel 

Each  sensor  is  equipped  with  5  IREDs 
covering  a  25°  range  of  deflection  data,  allow¬ 
ing  a  greater  resolution  of  data  collection.  A 
control  algorithm  is  required  in  order  to  obtain 
clean  and  stable  data  from  the  sensors  [32],  for 
this  reason,  an  embedded  microcontroller  is 
placed  inside  the  wheel  with  the  sensors.  The 
microcontroller,  particularly  the  Arduino 
Mega  2560  R3,  controls  which  sensors  to  use, 
collects  data  with  a  timestamp  to  be  saved  to 
an  on-board  memory  card,  and  transmits  the 
data  wirelessly  in  intervals.  Because  of  the 
control  algorithm  of  the  sensors  and  the  8 -bit 
Analog  to  Digital  Converter  (ADC)  of  the 
Arduino  Mega,  the  sample  frequency  is  lim¬ 
ited  to,  at  most,  100  Hz.  A  triple  axis  acceler¬ 
ometer  is  also  placed  inside  the  wheel,  primar¬ 
ily,  to  promote  efficient  data  collection  (by 
choosing  to  measure  from  sensors  only  at  the 
contact  patch)  but  can  also  be  used  to  compute 
the  test  tire’s  camber  and  toe  angles,  respec¬ 
tively,  as  well  as  the  angular  velocity. 
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3.2.4  WITS  -  Wireless  Capabilities 

Due  to  the  rotational  nature  of  the  wheel  and 
sealed  enclosure  between  the  tire  and  rim, 
making  the  WITS  system  wireless  reduces  the 
complexity  of  the  hardware  configuration. 
This  approach  also  lowers  the  overall  cost  of 
the  system  by  avoiding  the  use  of  special 
pneumatically  sealed  feed  through  valves, 
wheel  rim  machining,  and  slip  rings.  Yet  it  in 
turn,  raises  the  complexity  of  the  software 
running  the  system  and  forces  the  system  to 
be  power  efficient  for  long  term  running. 

To  fulfil  the  data  transmission  feature 
of  the  WITS  system,  the  XBee  Series  1  mod¬ 
ule,  which  utilizes  the  ZigBee/802.15.4  proto¬ 
col,  is  chosen  for  its  sufficiently-large  broad¬ 
casting  range  (30  m)  and  relativity  low  power 
requirements  (45  mA  at  3.3  V)[33].  Another 
advantage  is  that  interfacing  an  XBee  module 
with  the  Arduino  Mega  is  also  relativity 
straightforward  and  heavily  supported  in  the 
online  mechatronics  community. 

The  essential  idea  behind  using  the 
XBee  module  is  to  control  what  actions  the 
Arduino  Mega  will  follow  such  as  collecting 
data  for  a  test  run,  reading  and  extracting  data 
from  the  memory  card,  and  also  changing 
parameters  of  the  running  Arduino  Mega  pro¬ 
gram. 

3.2.5  WITS -Power 

Again,  because  of  the  wireless  nature 
and  low-budget  funding  of  the  WITS  system, 
a  battery  pack  approach  was  taken  to  supply 
power  to  the  WITS  system.  This  method  of 
power  delivery  was  chosen  for  its  ease  of  in¬ 
tegration  to  the  system  and  to  avoid  the  clutter 
and  cost  of  a  slip  rig,  which  would  have 
forced  the  project  to  use  a  feed  through  pneu¬ 
matically  sealed  valve. 

The  basic  idea  for  the  power  delivery  is 
to  have  a  battery  pack  separately  connected  to 
the  XBee  module  and  will  have  a  normally 
open  connection  to  the  Arduino  Mega  so  that 
at  the  user’s  request,  the  XBee  can  be  com¬ 
manded  to  allow  power  flow  to  the  Arduino 
Mega  and  sensors.  Power  conservation  relays 
substantially  on  the  XBee  modules  which  will 
be  in  a  “sleep”  mode  while  the  WITS  system 


is  not  in  use;  this  mode  of  operation  draws 
less  than  50  pA  of  current[33].  The  XBee  is 
configured  to  check  for  data  transmission 
every  hour,  if  it  finds  nothing,  it  will  return  to 
“sleep”,  if  it  receives  data,  it  will  turn  the  Ar¬ 
duino  Mega  on,  via  a  relay,  and  send  the  re¬ 
ceived  data  to  the  Arduino  Mega  so  that  it  can 
run  the  program  requested  by  the  user.  The 
power  switching  duty  of  the  XBee  is  accom¬ 
plished  by  its  simple  digital  input/out  and  ana¬ 
log  input  capabilities;  another  advantage  of 
using  the  XBee  module. 

As  for  the  batteries  themselves,  the 
WITS  system  will  be  powered  by  alkaline 
batteries.  Alkaline  batteries  were  chosen  to 
keep  the  battery  system  simple  and  avoid  the 
need  of  monitoring  the  charge/discharge  cycle 
that  common  polymer  lithium  ion  batteries 
require.  To  take  make  the  WITS  system  more 
green,  rechargeable  AA  batteries  can  be  used 
instead  of  the  current  AA  alkaline  batteries. 
Lastly,  to  minimize  weight  offset  caused  by 
the  battery  packs,  they  will  be  evenly  distrib¬ 
uted  around  the  wheel. 

3.2.6  WITS  -  Sensor  Validation 

To  ensure  the  accuracy  of  the  tire  de¬ 
flection  measured  by  the  WITS  system,  a 
validation  test  was  conducted  at  the  Virginia 
Tech  Performance  Engineering  Research  Lab 
(PERL).  The  validation  testing  consisted  of 
placing  a  wheel  fitted  with  sensors  on  the 
PERL’s  Tire-coupled  Quarter- Vehicle  test  rig 
and  comparing  the  sensors  measured  tire  de¬ 
flection  versus  the  deflection  actuated  by  the 
test  rig.  This  is  the  same  test  rig  used  in  the 
early  rig  surface  deflection  test  of  Section 
3.2.1. 

Several  tests  were  administered  to  see  how 
well  the  sensors  measured  different  tire  de¬ 
flections,  such  as  a  simple  step  input  to  a 
more  dynamic  sinusoidal  wave  fluctuation  of 
the  tire,  as  shown  in  Fig.  44. 
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time  (s) 

Fig.  44:  2.5cm  1/4  Hz  Sinusoidal  Input  (RED)  with  WITS 
measurement  (BLUE) 

Fig.  44  illustrates  that  the  WITS  re¬ 
sponse  is  a  fairly  accurate  measurement  of  the 
sinusoidal  measurement  of  tire  deflection.  The 
amplitude  difference  of  the  WITS  system  can 
be  attributed  to  the  small  tread  deflection  ex¬ 
perienced  as  the  rigid  hydraulic  actuation  ap¬ 
plies  the  load  to  the  wheel.  Also,  the  offset  in 
time  can  be  attributed  to  a  slight  difference  in 
timestamps  as  both  measurements  were  taken 
from  different  systems. 

Note  that  the  rig  is  also  capable  of 
measuring  the  vertical  force  applied  to  the 
wheel;  however,  this  data  was  useless  due  to 
the  air  leak  of  the  wheel.  The  air  leak  was  due 
to  early  design  efforts  of  the  WITS  system  to 
feed  the  wires  of  the  sensors  through  the  rim. 

Additional  validation  testing  for  the 
WITS  system  will  be  conducted  on  the  rig 
itself  using  a  Light  Detection  and  Ranging 
(LIDAR)  system.  The  fiber  optic  Doppler 
LIDAR  system  works  by  utilizing  scattered 
light  to  collect  data  on  a  number  of  parame¬ 
ters,  particularly  for  this  test  the  relative  ve¬ 
locity  of  the  wheel  angular  and  translational 
velocity.  Knowing  this  calculation  will  allow 
an  accurate  measure  of  the  slip,  thus  giving  a 
clear  measure  of  the  deformed  rolling  radius. 

3.2.7  Measure  Max  Sinkage 

Measuring  the  instantaneous  maximum 
sinkage  of  a  rolling  tire  in  soil  is  a  great  pa¬ 
rameter  of  interest  for  model  validation;  de¬ 
spite  the  difficult  nature  of  measuring  this 
parameter,  as  opposed  to  the  slightly  different 


rut  depth.  This  project  aims  to  measure  max 
sinkage  with  new  instrumentation.  The  max 
sinkage  can  be  calculated  with  the  use  of  the 
Terramechanics  rig  string  potentiometer 
(measure  vertical  wheel  displacement),  the 
WITS  system,  and  with  the  addition  of  an 
ultrasonic  sound  sensor  to  measure  the  soil 
height  in  the  path  of  a  rolling  test  tire.  Fig.  45 
illustrates  how  the  max  sinkage  is  calculated 
if  the  wheel  displacement  (h),  the  soil  height 
(h initial),  and  the  deformed  rolling  radius  (R^>_ 
formed )  are  known. 


Fig.  45:  Diagram  for  measuring  max  sinkage  on 
Terrmechanics  Rig 

Using  said  parameters  allows  the  calculation 
of  the  max  sinkage  as  shown  in  Eq.  (62). 

z0  ~  (jl  ^ deformed  )  —  ^ initial  (62) 

3.2.8  Normal  Load  Controller 

As  mentioned  in  Section  4.3  the  variation  in 
the  applied  wheel  normal  load  while  testing 
needs  to  be  addressed.  This  is  a  result  of  pas¬ 
sively  applying  a  load  to  the  test  wheel  via 
two  triple-bellow  air  springs,  shown  in  Fig. 
41.  Before  each  test  run,  air  is  applied  until 
the  desired  normal  load  is  achieved;  however, 
as  the  wheel  rolls  over  the  test  path,  the  wheel 
can  sink  or  rise  in  the  soil,  thus  changing  the 
normal  load.  To  counter  the  varying  normal 
load,  an  active  control  system  is  being  imple¬ 
mented  to  the  system. 

A  proportional,  integral,  derivative 
(PID)  control  algorithm  and  an  air  propor¬ 
tional  pneumatic  control  valve,  Enfield  Tech¬ 
nologies’  Mid,  will  be  employed  to  keep  ac¬ 
tive  control  of  the  wheel  normal  load.  The 
basic  idea  of  the  setup  is  to  command  the  Mid 
to  either  supply  10  bar  of  air  to  the  air  springs, 
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or  exhaust  the  air;  meanwhile,  using  the  wheel 
hub  sensor’s  vertical  force  readings  as  feed¬ 
back  to  the  controller.  Given  the  valve’s  abil¬ 
ity  to  supply  air  at  a  flow  rate  of  approxi¬ 
mately  2300  LPM  (at  10  bar)  in  less  than  5.6 
ms  response  time  to  a  change  in  signal 
input[34],  the  valve  is  expected  to  perform 
well  enough  to  maintain  a  constant  normal 
load  on  the  wheel. 

4.5  Soil  Work 

The  novel  aspect  of  the  project  is  model¬ 
ling  the  interaction  of  a  tire  over  a  deformable 
soil,  specific  to  this  study,  silty  sand;  there¬ 
fore,  efforts  to  characterize  the  soil  model 
with  realistic  parameters  are  of  great  im¬ 
portance.  As  demonstrated  in  the  Section  2 
and  Table  3,  multiple  soil  parameters  are  re¬ 
quired  to  characterize  the  tire-soil  interaction. 
Most  of  these  parameters  can  be  obtained  with 
the  bevameter  technique  [5],  while  others  re¬ 
quired  additional  experimental  testing  and  a 
drive  cylinder.  Given  the  difficulty  in  obtain¬ 
ing  a  machine  to  follow  the  bevameter  ap¬ 
proach  as  well  as  performing  the  various 
amounts  of  additional  experimental  soil  work, 
the  AYDL  has  decided  to  subcontract  Schna¬ 
bel  Engineering,  a  geotechnical  engineering 
company,  to  extract  soil  properties  of  our 
sandy  silt  soil. 

After  consulting  with  one  of  Schnabel’s 
associates,  C.J.  Smith,  about  the  intended  soil 
characterization  of  the  project,  the  level  of 
complexity  and  time  required  for  accurate  soil 
parameter  extraction  turned  out  to  be  much 
greater  than  expected.  Also,  the  parameters  of 
interest  to  the  geotechnical  world,  (commonly 
civil  engineers)  varied  greatly  from  the  pa¬ 
rameters  needed  from  a  terramechanics  stand¬ 
point.  For  this  reason  the  AVDL  decided  to 
launch  an  additional  project  to  work  alongside 
this  modelling  and  validation  project,  on  fur¬ 
ther  soil  parameter  extraction  and  transfor¬ 
mation  to  the  soil  model. 

Nonetheless,  lab  soil  testing  has  been  con¬ 
ducted  on  the  same  soil  and  soil  conditions 
that  the  Terramechanics  rig  will  keep  while 
performing  tire-soil  testing.  To  begin  this 


process,  Schnabel  had  to  measure  the  soil 
moisture  and  density  that  is  expected  for  ideal 
testing  conditions  via  their  nuclear  density 
gauge  [35]  in  order  to  replicate  the  same  sce¬ 
nario  in  the  laboratory  testing.  The  following 
subsection  will  explain  the  process  of  prepar¬ 
ing  the  soil  for  ideal  test  conditions  and  the 
subsequent  section  will  present  some  of  the 
laboratory  results  of  Schnabel. 

3.2.9  Soil  Preparation  and  Measurement 

To  ensure  repeatable  environmental  test 
conditions,  a  four  step  soil  preparation  proce¬ 
dure  was  established.  The  main  purpose  of  the 
soil  preparation  procedure  is  to  remove  any 
standing  history  of  soil  deformity,  create  a 
uniform  level  of  soil,  and  compact  the  soil  to  a 
desired  density.  Multiple  soil  preparation  pro¬ 
cedures  have  been  developed  over  the  years, 
such  as  the  automated  “Mixing-and- 
Compaction”  device  in  [36]  and  semi- 
automated  procedure  of  [20].  Investigations  of 
the  impact  the  soil  reconditioning  procedure 
itself  would  have  on  the  study  results  have 
been  explored  as  well  in  [37,  38],  The  proce¬ 
dure  followed  in  the  AVDL  is  similar  to  that 
of  [20]  but  has  been  modified  to  cater  to  the 
different  soil  structure  of  silty  sand.  The  pro¬ 
cess  begins  with  tilling  the  soil  using  an  elec¬ 
trical  tilling  machine;  three  passes  along  a  3  ft 
stretch  are  prepared  to  completely  loosen  the 
soil.  The  soil  surface  is  then  raked  to  provide 
a  uniform  and  relativity  smooth  surface.  Next 
the  Terramechanics  rig  carriage  is  fitted  with 
a  smoothing  blade  and  thoroughly  flattens  out 
the  soil  surface  to  a  level,  smooth  surface. 
The  final  step  requires  making  two  passes 
with  a  variable-weight  lawn  roller  to  compact 
the  soil. 

Following  the  soil  preparation  procedure, 
the  soil  compaction  and  moisture  is  measured. 
A  hand-held  cone  penetrometer  is  the  used  to 
measure  the  compaction  versus  sinkage.  The 
moisture  content  is  measure  using  a  simple 
ASTM  standard  [39]  of  microwaving  the  soil 
to  measure  the  weight  difference  attributed  to 
evaporated  water.  Averages  of  these  meas¬ 
urements  are  taken  to  confirm  overall  con- 
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sistent  soil  parameters  along  the  entire  test 
path  of  prepared  soil. 

3.2.10  Soil  Parameter  Extraction 

Under  the  advisory  of  C.J.  Smith  and 
Jody  Priddy,  a  Ph.  D  student  in  the  Geotech¬ 
nical  Engineering  Program  of  the  Civil  &  En¬ 
vironmental  Engineering  Dept,  of  VT,  the  test 
methods  that  would  benefit  the  characteriza¬ 
tion  of  soil  in  the  tire  model  were  discussed 
and  executed.  To  begin  the  laboratory  testing, 
it  was  important  to  understand  the  physical 
makeup  of  the  soil,  which  was  initially  be¬ 
lieved  to  be  “sandy  loam”  but  turned  out  to  be 
more  of  silty  sand,  according  to  the  civil  engi¬ 
neering  and  highway  construction  classifica¬ 
tion  system  [40,  41].  The  following  test  per¬ 
formed  was  the  standard  proctor  method 
which  is  used  to  understand  the  moisture- 
density  relation  of  the  soil  [42],  This  test  re¬ 
vealed  the  maximum  point  of  soil  density  ver¬ 
sus  moisture  content.  Similar  to  a  cone  pene¬ 
trometer,  a  California  Bearing  Ratio  (CBR) 
test  was  administered  to  measure  the  relative 
resistance  to  penetration  of  the  soil  [43],  The 
results  of  this  data  can  be  closely  correlated 
with  the  cone  penetration  data  obtained  during 
testing  in  the  Terramechanics  rig.To  investi¬ 
gate  the  shear  strength  of  the  soil,  a  triaxial 
compression  test  was  implemented.  The  test 
was  setup  to  measure  the  soil’s  unconsolidat¬ 
ed,  undrained  (UU)  internal  shear  strength 
[44],  again,  to  mimic  the  soil  test  conditions 
in  tire-soil  testing.  The  test  was  administered 
on  three  different  specimens  of  varying  con¬ 
fining  stress  (3  kPa,  10  kPa,  and  21  kPa). 
Among  the  resulting  data  of  triaxial  testing  of 
the  three  specimens  the  obtained  soil  stresses 
were  plotted  in  Mohr  circles  in  shear  versus 
normal  stress  in  the  failure  stress  states;  out¬ 
lining  a  Mohr-Coulomb  failure  envelope  de¬ 
fined  by  an  internal  friction  angle  of  17  de¬ 
grees  and  cohesion  of  9.7  kPa.  Fig.  46 
demonstrates  the  resulting  Mohr  circle  plot. 


0  0  14  28  41  55  69 


Total  Normal  Stress  (kPa) 

Fig.  46:  Mohr  Stress  Circles  at  2.5%  Axial  Strain  Criterion; 
from  Schnabel  lab  report 

To  measure  the  compressibility  behaviour  of 
the  soil,  a  one-dimensional  consolidation  test 
was  performed  [45],  In  this  test,  the  soil  was 
subjected  to  incremental,  monotonic  static 
from  12  kPa  up  to  a  maximum  of  766  kPa, 
with  unload-reload  hysteresis  cycles  at  192 
kPa  and  766  kPa. 

5  Conclusions 

This  paper  presents  a  proposed  semi- 
analytical  soft  soil  tire  model  and  its  experi¬ 
mental  validation  methodology.  The  model  is 
intended  for  traction,  ride,  handling,  and  du¬ 
rability  applications.  The  work  on  the  project 
is  ongoing  and  simulation  results  as  well  as 
validation  tests  will  be  presented  in  future 
publications. 
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7  Nomenclature 


sd  Longitudinal  slip 


b  Tire  width 


Vx  Vehicle  longitudinal  velocity 


R  Initial  radius  approximation  di- 

R  mensionless  variable 

c  Soil  cohesion 

Experimental  constant  used  to 
C°  determine  0m 

Ci  Multipass  experimental  parameter 

csd  Multipass  experimental  parameter 

Experimental  constant  used  to 
Cm  determine  0m 

8  Tire  deformation 

,  Longitudinal  shear  deformation 

x  modulus 

ky  Lateral  shear  deformation  modulus 

kxi  Multipass  experimental  parameter 


ksd 

k0 

fci 

k2 

h 

k4 

^5 

K 

Cl 


c2 


c4 

Ri 


Multipass  experimental  parameter 

Static  modulus  of  soil  deformation 

Sidewall  radial  spring  stiffness  (in¬ 
plane) 

Wheel-sidewall  radial  spring  stiff¬ 
ness  (in-plane) 

Inter-element  radial  spring  stiff¬ 
ness  (in-plane) 

Tangential  inter-element  spring 
stiffness  (in-plane) 

Lateral  inter-element  spring  stiff¬ 
ness  (out-of-plane) 

Radial  inter-element  spring  stiff¬ 
ness  (out-of-plane) 

Sidewall  radial  damping  (in-plane) 

Sidewall  circumferential  damping 
(in-plane) 

Lateral  inter-element  damping 
(out-of-plane) 

Wheel-sidewall  radial  damping 
(in-plane) 

Loaded  radius 
Un-deformed  radius 


Rdeformed  Deformed  rolling  radius 


z0  Maximum  sinkage 

n  Pressure-sinkage  index 

m  Penetration  velocity  exponent 

mbead  Bead  mass 

mi  bt  Belt  plus  tread  element  mass 

^thread  Tread  mass 

mbeit  Bed  mass 

msidewaii  Total  sidewall  mass 
mi  s  Sidewall  element  mass 

mwheei  Wheel  mass 

Nm  Number  of  masses 

an  Normal  stress 

tx  Longitudinal  shear  stress 

r y  Lateral  shear  stress 

jx  Longitudinal  shear  displacement 

jy  Lateral  shear  displacement 

p  Radial  pressure 

pt  Tire  internal  pressure 

0  Soil  angle  of  internal  friction 

y  Soil  density 

Ys  Unit  weight  of  soil 

Ysi  Multipass  experimental  parameter 

Ysd  Multipass  experimental  parameter 

Ny  Soil  specific  weight  coefficient 

Nc  Soil  cohesion  coefficient 

Nq  Soil  surcharge  load  coefficient 

Surcharge  load  from  accumulated 
^  bulldozed  soil 

z  Sinkage 

7  Initial  radius  approximation  di- 

mensionless  variable 
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z  Vertical  penetration  velocity 

z0  Maximum  sinkage 

~  Angle  of  friction  between  soil  and 

Of 

1  tire 

^  Central  angle  describing  mass  po¬ 

sition 

6a  Angle  of  transition 

0b  Trailing  Edge  Angle 

9e  Entry  angle 

df  Leading  edge  angle 
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9  Appendix 

A.  Forces  acting  on  wheel  due  to  springs  and  dampers 

The  forces  acting  on  the  wheel  due  to  the  springs  and  dampers: 
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where  the  first  index  denotes  the  direction  of  the  force,  and  the  second  one  the  mass  number.  Note 
that  the  mass  number  index  0  refers  to  the  rigid  wheel. 


B.  Equations  of  motion  on  the  Side  wall 

The  forces  acting  on  the  springs  and  dampers  connected  to  the  left  sidewall  are  presented  in  equa¬ 
tion  (9-2),  (9-3)  and  (9-4). 
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F4fu  =  fc4(— Zj-u  sin  Ad  +  Xi_lti  cos  Ad  —  x^)  —  k4(— zi+1>l  sin  Ad  —  xi+1>[  cos  Ad  +  Xif) 

+  c4{-Zi_u  sin  Ad  +  Xi-U  cos  Ad  -  xM)  (9-3) 

—  c4(— zi+1i  sin  Ad  —  xi+1>i  cos  Ad  +  xtf) 


F3  i  i  =  k3 {xi_1:i  sin  Ad  +  Zi_lit  cos  Ad  —  zi  t )  —  k3 ( xi+l  i  sin  Ad  —  zi+l  i  cos  Ad  +  ztf) 


(9-4) 


On  the  other  hand,  the  forces  on  the  right  sidewall  are  the  following, 
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F4ir  =  kA{—Zi_l  r  sin  Ad  +  Xi_lr  cos  Ad  —  xiT) 

—  k4(—zi+lr  sin  Ad  —  xi+l  r  cos  Ad  +  xi  r)  (9_6) 

+  c4(— Zj_l  r  sin  Ad  +  ij_ l  r  cos  Ad  —  xi  r) 

—  c4(— zi+1  r  sin  Ad  —  xi+l  r  cos  Ad  +  xi  r) 

F3  i  r  =  k3 (xj_l  r  sin  Ad  +  Zj_x  r  cos  Ad  —  zi  r)  —  k3  ( xi+l  r  sin  Ad  —  zi+l  r  cos  Ad  +  zi  r)  (9_7) 


C.  Equation  of  motion  Belt  and  Tread 

The  forces  acting  on  the  belt  and  tread  elements  due  to  the  connected  springs  and  dampers  are  de¬ 
scribed  below. 
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As  such,  the  equations  of  motion  for  the  belt  and  tread  layer  are  defined  by  equation 
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D.  Shear  Stress 


1  —  e 


Rj[(0g-0e)-(cos  fle- cos  6a)+Srf(cos  ee-cos  fla)]  \  2  /  ftt(l~Sd)(0q-0e)  tan  «c 


+  1  —  e 


=  1 


(9-10) 


Tx.iSiWi 


ee>  9  >ea 


Fxcp,i  ~  ) 


FxSd 


2 

(nytanac)  +  (pxsd)2 


T-y,isiwi 


9a~>  9  >  6r 


9e>  9  >9a 


Fycp.i  ~  A 


liy  tan  ac 


2 

(jUytanac)  +  (p.xsdy 


:Anax,tS;Wj  9a  >  9  >  9r 


(9-11) 


(9-12) 


where  Ux  and  fly  are  the  friction  coefficients;  and  they  are  computed  using  Eq.  (48). 
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